1/2 


AD-A122  127  INTERACTIVE  MECHANISMS  OF  SLIDING-SURFACE  BEARINGS(U) 

DAVID  W  TAVLOR  NAVAL  SHIP  RESEARCH  AND  DEVELOPMENT 
CENTER  BETHESDA  MD  T  L  DAUGHERTV  ET  RL.  AUG  82 
UNCLASSIFIED  DTNSRDC-82/119  F/G  11/8  NL 


OF  SLIDING-SURFACE  BEARINGS  DTNSR DC-82/1  IS 


DAVID  W.  TAYLOR  NAVAL  SHIP 
RESEARCH  AND  DEVELOPMENT  CENTER 

litlwidi,  Maryland  20084 


INTERACTIVE  MECHANISMS  OF 
SLIDING-SURFACE  BEARINGS 


Thomas  L.  Daugherty,  DTNSRDC 


Coda  H.  T.  Pan,  Columbia  University 


APPROVED  FOR  PUBLIC  RELEASE: 
DISTRIBUTION  UNLIMITED 


SHIP  MATERIALS  ENGINEERING  DEPARTMENT  1J  \ 
RESEARCH  AND  DEVELOPMENT  REPORT 


0  \983 


DTNSR  DC-82/1  IB 


m 


,83  C9  28  09  6 


f- 


MAJOR  DTNSRDC  ORGANIZATIONAL  COMPONENTS 


UNCLASSIFIED 


SECURITY  CLASSIFICATION  OF  THIS  PAGE  (When  Data  Enter ed) 


REPORT  DOCUMENTATION  PAGE 


report  number 
DTNSRDC-82/119 


4.  TITLE  (end  Subtitle) 


12.  GOVT  AC 


INTERACTIVE  MECHANISMS  OF  SLIDING 
SURFACE  BEARINGS 


READ  INSTRUCTIONS 
BEFORE  COMPLETING  FORM 


RECIPIENT'S  CATALOG  NUMBER 


S.  TYPE  OF  REPORT  &  PERIOD  COVERED 


Research  and  Development 


6.  PERFORMING  ORG.  REPORT  NUMBER 


7.  AUTHOR/*.) 

Thomas  L.  Daugherty,  DTNSRDC 

Dr.  Coda  H.  T.  Pan,  Columbia  University 


9.  PERFORMING  ORGANIZATION  NAME  AND  AODRESS 

David  Taylor  Naval  Ship  R&D  Center 
Annapolis,  MD.  21402 


II.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

Naval  Sea  Systems  Command  (SEA  05R25) 
Washington,  D.  C.  20362 


4.  MONITORING  AGENCY  NAME  A  AODRESS/lf  different  from  Controlling  Olllco) 


•  ■  CONTRACT  OR  GRANT  NUMBER/*) 


10.  PROGRAM  ELEMENT.  PROJECT,  TASK 
AREA  A  WORK  UNIT  NUMBERS 

Program  Element  62761N 
Task  Area  SF  61-541-502 


12.  REPORT  DATE 

August  1983 


IS.  NUMBER  OF  PAGES 

141 


IS.  SECURITY  CLASS,  (ol  thle  report) 


16.  DISTRIBUTION  STATEMENT  (ol  thle  Report) 


UNCLASSIFIED 


15a.  DECLASSIFICATION/ DOWNGRADING 
SCHEDULE 


APPROVED  FOR  PUBLIC  RELEASE:  DISTRIBUTION  UNLIMITED. 


17.  DISTRIBUTION  STATEMENT  (of  the  ebetrmct  entered  In  Block  20,  H  different  from  Report) 


19.  KEY  WORDS  (Continue  on  teveree  aide  If  neceeeery  end  Identity  by  block  number) 

Hydrodynamic  Surface  Roughness 

Mixed  Lubrication  Friction 

Fluid  Film  Bearings 
Cavitation 


20.  ABSTRACT  (Qbntlnue  on  reeeree  elde  If  neceeeery  and  Identity  by  block  number) 

_ The  theory  of  hydrodynamic  lubrication^as  been  adapted  in  this  report 

<to  include  surface  roughness  effects  and  asperity  interactions.  This  theory 
becomes  a  means  of  studying  the  behavior  of  sliding-surface  bearings  in 
transition  from  hydrodynamic  into  mixed  lubrication.  A  model  problem,  in  the 
form  of  a  crowned  tilt-pad  bearing  with  longitudinal  roughness,  was  solved. 
In^the  treatment  of/the  model  problem,  end-leakage  was  neglected  and  a 
constant  lubricant  viscosity  was  assumed.  Use  of  numerical  computational 


I  JAN  71  14/3  Y  COITION  OF  I  NOV  6*  II  OBtOLCTC 
S/N  0102.014*6601 


UNCLASSIFIED 

SECURITY  CLASSIFICATION  OF  TNIS  PAOC  (When  Date  Entered) 


UNCLASSIFIED 


-LLUHITY  CLASSIFICATION  OF  THIS  PAGErOTian  Data  En tired) 


(Block  20  continued) 

y 

,y  ' 

-^►methods  could  readily  remove  these  restrictions. 

In  real-CllW^operating  environments  of  ship  components  and  shipboard 
machinery,  the  lubrication  process  is  made  more  complex  by  thermo- 
elastohydrodynamic  effects  and  the  likely  occurrence  of  film  striation. 
Because  these  complicated  phenomena  are  only  qualitatively  understood,  no 
adequate  or  comprehensive  method  exists  for  designing  critical  machine 
components  such  as  face  seals,  elastomeric  sterntube  bearings,  rudder  stock, 
diving  plane  bearings,  and  main  shaft  thrust  bearings.  'r‘"e.  t'p--  T  . 

The  work  presented  herein -is  part  of  a  Wrt  comprehensive  plan  to 
develop  a  better  understanding  of  the  interactive  wear  mechanisms  of 
sliding  surfaces.  The  theoretical  basis  is  the  application  of  Reynolds' 
equation  to  the  mixed  lubrication  regime.  A  closed-form  solution  of 
Reynolds'  equation  applied  to  a  crowned  tilt-pad  thrust  bearing  in  trans¬ 
ition  from  the  hydrodynamic  to  the  mixed  lubrication  region  is  presented 
to  demonstrate  some  of  the  principles  of  this  approach.^  Associated  topics 
,  related  to  the  extension  of  such  analyses  are  alsQ/discussed.  Concepts 
ff'-5,  developed  hum  in  are  applicable  to  a  wide  variety  of  bearings  and  seals.  > 

'  They  provide  a  basis  for  future  work  in  tribology  of  improved  machinery  ^ 
design  methods. 


Accession  For 

NTIS  GRAM 

DTIC  TAB 

Unannounced 

□ 

Justification _ 

By 

-  ■  — 

Distribution/ 
Availability  Codes 

(Avail  and/or 
Spec ial 


UNCLASSIFIED 


SECURITY  CLASSIFICATION  OF  THIS  F AOE(HTi»n  Data  Entered) 


TABLE  OF  CONTENTS 


Page 


LIST  OF  FIGURES .  iv 

LIST  OF  TABLES .  iv 

NOMENCLATURE  .  v 

ABSTRACT  .  1 

ADMINISTRATIVE  INFORMATION  .  1 

INTRODUCTION  .  2 

BACKGROUND  .  4 

THEORETICAL  TOPICS  RELATED  TO  THE  STUDY  OF  MIXED  LUBRICATION  OF 

SLIDING  SURFACE  BEARINGS .  10 

1.  FILM  STRIATION .  10 

Incipience  Point  of  Film  Rupture  .  10 

Viscous  Shear  of  Striated  Film .  12 

2.  THERMO-ELASTOHYDRODYNAMIC  EFFECTS  .  13 

General  Discussion .  13 

Conservation  Laws  for  Fluid  Film .  16 

Thermal  Effects . 26 

Analysis  of  the  Temperature  Field  .  28 

3.  COMPUTATION  TECHNIQUES  .  35 

Finite  Difference  Computation  of  the  Infinitely  Long  Slider  .  36 

Description  of  Thrust  Bearing  Film  Profile.  ...  .  42 

High  Accuracy  Numerical  Algorithm  for  Computing  Fluid 

Film  Pressure  .  .......  .  44 

4.  SURFACE  ROUGHNESS  EFFECTS  ...  .  55 

Mixed  Lubrication  of  Sliding  Surface  Bearings  with 

Longitudinal  Roughness .  57 

Stochastic  Computations  for  Mixed  Lubrication  with  Skewed 

Roughness  Distribution  Function .  70 

CONCLUSIONS .  78 

RECOMMENDATIONS  .  78 


ACKNOWLEDGEMENTS 


80 


REFERENCES 


APPENDIX  A  -  THE  CROWNED,  MOMENT  BALANCED  TILT  PAD  BEARING 

IN  MIXED  LUBRICATION . 

APPENDIX  B  -  INTEGRALS  FOR  CIRCUMFERENTIAL  AND  AXIAL  COMPUTATION 

APPENDIX  C  -  INTEGRALS  FOR  RADIAL  COMPUTATION  . 


LIST  OF  FIGURES 


1  -  Stribeck  Curve  Showing  the  Three  Main  Lubrication  Regimes  . 

2  -  Two-Dimensional  Tilt  Pad  Bearing  .  . 

3  -  Christensen's  Longitudinal,  One-Dimensional  Roughness  Model 

4  -  Crowned  Tilt-Pad  Bearing  with  Longitudinal,  One-Dimensional 

Roughness . .  . 

5  -  Thin  Film  Between  Nearly  Parallel  Surfaces  . 

6  -  Schematic  Thermal  System  . 

7  -  Sector  Pad  of  a  Thrust  Bearing  ..  . 

8  -  Schematic  of  a  Nonuniform  Mesh  Cluster  . 

9  -  Applicable  Domain  for  Ensemble  Averaging  . 


LIST  OF  TABLES 

1  -  Thermo-Elastohydrodynamic  Phenomena  in  Typical  Naval 

Machine  Elements  . 

2  -  Minimum  Recommend  Mesh  Points  . 


3  -  Dependence  of  Third  Moment  and  Square  of  Half  Maximum 
Height  on  the  Skew-Bias  Parameters  ...  .  , 


NOMENCLATURE 


a0»  al*  a2 
Ac 

A,  B,  C,  D 
Ai>  Bt 

Ar»  Az 

[*<0)1,,  [*<->],. 

bO»  bl»  b2 
[B,] 

b 


cO»  Cl»  C2 

c 

CE 


D/Dt 

Ee 

E,  F 
E  {  } 
Eo  {  } 
AE  {  } 


Coefficients  in  polynominal  expansion  of  E{h3} 

Area  fraction  of  asperity  contact 
Corners  of  sector  pad 
Coefficients  at  ith  mesh  point 
Arbitrary  vector 

(Radial,  circumferential t  axial)  component  of  A 
Matrices  in  conduction  problem  [equation  (2-34)| 
Coefficients  in  polynominal  expansion  of  6qUE{H} 
Influence  coefficient  matrix  of  {q(}  on  {T}i 
Width  of  bearing  pad 

Half  total  range  of  the  random  film  thickness 
variable,  8 

c/o 

Coefficients  of  polynomial  expansion 
Nominal  bearing  clearance 
Constant  of  integration 

Specific  heat  at  a  constant  volume  process 
Coefficients  at  ith  mesh  point 
Substantive  time  derivative 
Specific  internal  energy 
Mid-points  of  radial  sector  sides 
Expectation  or  stochastic  average  of  {  } 
Uninterfered  expectation  of  (  ) 

Expectation  deficiency  of  (  )  due  to  interference 


Frictional  force 


Normal  force 


t  mX  ’j  TSTTi  T*  ”•  v  ■> 


Fn 

f 

F  (x,  y/h) 
G 

M 

h 

hj  ,  h2 

ha 

hb 

he 

ha 

he 

hf 

h3 


fy 


Jc 

!p 

In 


A  stochastic  function 

Temperature  profile  as  a  function  of  x  and  y/h 

Symmetrical  polynomial  generating  function  for 
stochastic  computations 

Film  thickness,  including  roughness  deviations 

Smooth  fluid  film  thickness,  h  =  h2  -  hj 

Distance  to  (lower,  upper)  bearing  surface  from 
reference  bearing  surface 

Circumferential  crown 

Radial  crown 

Crown  increment 

Edge  pitching 

Edge  rolling 

Radial  helix 

Longitudinal  roughness  deviation  function  of  one 
surface 

Heat  transfer  coefficient  at  exposed  surface 

Mesh  point  index,  ranges  from  1  to  N+l 

Unit  base  vectors  of  the  cylindrical  polar  coordinate 
system 

kth  moment  integral  of  inverse  viscosity;  k  may  be  0, 
1 ,  or  2 

Viscosity  weighted,  Couette  velocity  profile  function 

Viscosity  weighted,  Poiseuille  velocity  profile 
function 


Jc»  Jp 


nth  moment  integral  of  G 
Integral  of  (Ic,  Ip) 


Viscosity  weighted  (Couette,  Poiseuille)  temperature 
transport  integral 

Bearing  pad  length 

Coefficents  in  Walther-ASTM  formula  for  temperature 
effect  on  viscosity 

Mid-points  of  peripheral  sector  arcs 

Unit  vector  perpendicular  to  the  nominal  bearing 
surface 

Total  number  of  mesh  intervals 
Origin 

Pressure  of  fluid  film 
Cavitation  pressure 
Pressure  of  fluid  film  at  xj 
Center  point  of  sector 

Histogram  of  probability  density  function  of  roughne 
deviations  of  (lower,  upper)  bearing  surface 

Probability  distribution  function  of  hs 

x  -  component  of  flow  flux 

z  -  component  of  flow  flux 

Conduction  heat  flux 

Column  vector  representing  conduction  heat  flux  at 
bearing  surface  by  body  I 

Specific  internal  dissipation 

Radial  coordinate  in  polar  notation 

Inner  radius  of  sector 

Mean  radius  of  sector 

Outer  radius  of  sector 


Nominal  radius 


Exposed  surface  of  bearing 
Time 


Fluid  temperature 
Ambient  temperature 

Column  vector  representing  temperature  in  body  I  at  yj 
Temperature  of  exposed  surface 
x  Component  of  fluid  velocity 
Sliding  speed 

Fluid  velocity  component  tangential  to  the  nominal 
bearing  surface 

x  component  of  U  at  (lower,  upper)  bearing  surface 

Fluid  velocity  vector 

— ► 

Reference  value  of  V 
Thickness  of  body  I 

z  component  of  U  at  (lower,  upper)  bearing  surface 

Coordinate  of  the  cavitation  or  rupture  boundary 

Coordinate  of  ith  mesh  point 

Cartesian  coordinates 

Discretized  y  coordinate 

Dummy  variable  of  integration  for  y 

Arc  angle  of  sector 

Bulk  coefficient  of  thermal  expansion 
Roughness  deviation 
Skew-bias  of  roughness  deviation 
8  o/o 


vi  i  i 


slts2 

8  (  ) 

Ax 

e 

e 

0s 

0‘s 

K 

Kl 

M 

Mo 

V 

P 

a 

Oj  ,  o2 

i 

<t> 

4> 

4*o  *  4*j  »  4*2 

V 

(i ) 

"l»  w2 


V 


Roughness  deviations  of  (lower,  upper)  bearing 
surface,  (hj-Hj,  h2~H2) 

Increment  of  (  ) 

Mesh  interval 

Film  integral  of  temperature 

Circumferential  angular  coordinate 

Root  mean  third  moment  stochastic  integral  of  8 

9s/o 

Thermal  conductivity  of  fluid 
Thermal  conductivity  of  body  I 
Viscosity  coefficient  of  fluid 
Reference  viscosity  coefficient 
Kinematic  viscosity  coefficient  of  fluid 
x  -  X; 

Fluid  density 
Combined  roughness 
Rms  value  of  (8j,  82) 

Viscous  shear  stress 

Leakage  flux  through  bearing  surfaces 

Film  flux  (expectation  of) 

Coefficients  in  polynomial  expansion  of  4> 

Film  flux  vector 
Angular  rotational  speed 

Angular  rotational  speed  of  (lower,  upper)  bearing 
surface 

Three-dimensional  gradient  operator 


A 


IX 


ABSTRACT 


The  theory  of  hydrodynamic  lubrication  has  been  adapted 
in  this  report  to  include  surface  roughness  effects  and 
asperity  interactions.  This  theory  becomes  a  means  of 
studying  the  behavior  of  sliding-surface  bearings  in 
transition  from  hydrodynamic  into  mixed  lubrication.  A  model 
problem,  in  the  form  of  a  crowned  tilt-pad  bearing  with 
longitudinal  roughness ,  was  solved.  In  the  treatment  of  the 
model  problem,  end-;  'akage  was  neglected  and  a  constant 
lubricant  viscosity  was  assumed.  Use  of  numerical 
computational  methods  could  readily  remove  these 
restrictions . 

In  real-life  operating  environments  of  ship  components 
and  shipboard  machinery,  the  lubrication  process  is  made  more 
complex  by  thermo-elastohydrodynamic  effects  and  the  likely 
occurrence  of  film  striation.  Because  these  complicated 
phenomena  are  only  qualitatively  understood,  no  adequate  or 
comprehensive  method  exists  for  designing  critical  machine 
components  such  as  face  seals,  elastomeric  sterntube 
bearings,  rudder  stock,  diving  plane  bearings,  and  main  shaft 
thrust  bearings. 

The  work  presented  herein  is  part  of  a  more 
comprehensive  plan  to  develop  a  better  understanding  of  the 
interactive  wear  mechanisms  of  sliding  surfaces.  The 
theoretical  basis  is  the  application  of  Reynolds'  equation  to 
the  mixed  lubrication  regime.  A  closed-form  solution  of 
Reynolds'  equation  applied  to  a  crowned  tilt-pad  thrust 
bearing  in  transition  from  the  hydrodynamic  to  the  mixed 
lubrication  region  is  presented  to  demonstrate  some  of  the 
principles  of  this  approach.  Associated  topics  related  to 
the  extension  of  such  analyses  are  also  discussed.  Concepts 
developed  herein  are  applicable  to  a  wide  variety  of  bearings 
and  seals.  They  provide  a  basis  for  future  work  in  tribology 
of  improved  machinery  design  methods. 

ADMINISTRATIVE  INFORMATION 

This  report  covers  work  conducted  under  the  Marine  Tribology  Block  Program 
(PE62761N,  Task  Area  SF  61-541-502,  Work  Unit  2832-101).  The  coauthor,  Dr.  Coda 
H.  T.  Pan,  is  a  consultant  with  DTNSRDC  under  Contract  N00167-81-0149.  Dr.  Pan 
a  Professor  of  Mechanical  Engineering  at  Columbia  University.  The  Marine 
Tribology  Block  program  is  sponsored  by  Dr.  H.  H.  Vanderveldt  of  NAVSEA  05R25. 
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INTRODUCTION 


Most  sLiding-surface  bearings  are  designed  to  operate  with  the  bearing  fully 
separated  from  its  mating  surface  by  a  lubricant  film.  These  bearings  are 
referred  to  as  self-acting  or  hydrodynamic  bearings.  The  lubricant  film  develops 
load  capacity,  primarily  due  to  a  wedge  effect  between  the  moving  surfaces. 

Various  lubrication  modes  of  such  bearings  can  be  explained  by  considering  the 

Stribeck  curve  shown  in  Figure  1.  The  Stribeck  curve  represents  the  general 

1  ☆  » 

characteristics  of  lubricated  moving  surfaces  as  a  function  of  the  lubricant 
viscosity,  p,  the  velocity,  V,  and  the  normal  load,  Fn  (or  pressure,  p).  As 
Figure  1  shows,  there  are  three  main  lubrication  regimes:  hydrodynamic,  mixed, 
and  boundary.  The  regime  of  operation  depends  upon  the  geometry,  the  materials, 
the  operating  conditions,  and  the  relative  geometry  of  the  opposing  surfaces. 

In  the  hydrodynamic  regime  a  wedge  of  lubricant  provides  a  thick  enough 
lubricant  film  to  keep  the  surfaces  completely  separate,  such  that  their  combined 
roughness,  o,  is  substantially  smaller  than  the  local  film  thickness,  h.  The 
frictional  resistance  results  from  viscous  shearing  of  the  lubricant.  Since  the 
surfaces  are  separated,  no  contact  wear  takes  place.  Under  these  conditions, 
analytical  techniques  for  the  study  of  bearing  performance  are  well  developed. 

Mixed  lubrication  begins  when  the  local  film  thickness,  h,  approaches  that  of 

the  combined  surface  roughness,  o,  of  the  bearing  and  the  mating  surfaces.  The 

pressure  generated  by  the  hydrodynamic  wedge  effect  and  by  the  asperity  contact 
between  the  bearing  and  the  mating  surface  share  the  load  support.  Wear  occurs  at 

the  contact  points.  Frictional  resistance  is  a  result  of  these  asperity 

interactions  plus  shearing  of  the  lubricant. 

As  load  is  increased  and/or  speed  or  lubricant  viscosity  are  decreased,  the 
coefficient  of  friction  and  the  wear  rate  increase  rapidly.  Eventually  all  load 
will  be  supported  by  asperity  contacts.  This  is  the  boundary  lubrication  region. 
Physico-chemical  interactions  at  the  solid-lubricant-solid  interface  determine  the 
friction  and  wear  behavior. 

In  a  wide  variety  of  sliding  surface  bearing  and  seal  applications  in  ships, 
not  enough  lubricant  film  exists  under  certain  operating  conditions  to  completely 


*Ref erences  are  listed  on  page  81. 
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Figure  1  -  Stribeck  Curve  Showing  Che  Three  Main  Lubrication  Regimes.  The 
Characteristics  of  Moving  Lubricated  Surfaces  are  a  Function  of  Lubricant 
Viscosity,  n,  the  Velocity,  V,  and  the  Normal  Load,  FN;  h  is  the  Fluid  Film 
Thickness  and  o  the  Combined  Roughness  of  the  Surfaces. 
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separate  the  mating  surfaces.  Such  conditions  cause  mechanical  wear.  Examples 
are  ship  main  thrust  bearings,  lineshaft  journal  bearings,  and  compliant-surfaced 
sterntube  and  strut  bearings.  For  these  bearings,  films  are  often  inadequate  at 
slow  shaft  speeds  and  under  transient  conditions  of  starts  and  stops.  Face  seals 
are  unique  in  having  to  operate  most  of  their  lives  in  a  contact  and  wear  mode, 
since  they  are  designed  to  minimize  leakage.  A  good  face  seal  design  thus 
successfully  controls  wear  during  normal  operation  without  adequate  fluid  film. 

Even  though  the  design  methodology  for  full-film  sliding-surface  bearings  is 
relatively  well  established,  tribology  laboratories  regularly  perform  friction  and 
wear  tests  in  search  of  better  bearing  materials.  The  challenge  is  to  improve 
service  reliability  of  ever  more  complex  mechanical  systems.  A  comprehensive 
knowledge  of  the  interactive  wear  mechanisms  of  sliding  surfaces  is  of  paramount 
importance.  Because  of  the  multitude  of  phenomena  involved,  an  interdisciplinary 
point  of  view  must  be  maintained.  Friction  and  wear  studies,  which  are 
traditionally  pursued  with  empirical  methods,  can  and  must  be  vitalized  by 
computer-aided  techniques  of  design,  analysis,  and  research.  Theory  for  mixed 
lubrication  plays  a  pivotal  role  in  the  meaningful  study  of  interactive  wear 
mechanisms.  It  should  address  four  main  issues. 

o  How  can  one  predict  the  start  of  asperity  interactions? 

o  What  are  the  important  thermo-mechanical  consequences  of  asperity 
interactions? 

o  What  are  the  relevant  roles  of  lubricant  chemistry,  surface  chemistry,  and 
near-surface  behavior  of  materials? 

o  How  does  the  wear  process  alter  the  topographical  and  chemical  properties 
of  the  surface? 

Progress  in  finding  answers  to  these  questions  will  advance  design  technology  for 
machinery  systems  in  ships. 

BACKGROUND 

Reynolds'  equation  has  been  used  to  study  hydrodynamics  of  lubricating 
films  since  its  development  by  Professor  Osborne  Reynolds2  in  1886.  It  is  a 
differential  equation  for  the  pressure  distribution  in  a  lubricant  film  and  takes 
into  account  the  density  and  viscosity  of  the  lubricant,  the  velocities  of  the 


solid  surfaces,  and  the  shape  of  the  lubricant  film.  The  most  general  Cartesian 
form  of  Reynolds'  equation  can  be  written  as: 

i(iS£t)+i(*i).,|1i( ,*)♦*! P*( +  »  pA(w1+«r2) 

dx  \  n  dx  )  dz\pdzJ  l  at  V  /  v  1  2/  dz  \  1  2/1 


where:  h 


=  smooth  fluid  film  thickness 


p  =  pressure  of  fluid  film 
x,z  =  cartesian  coordinates 

p  =  fluid  density 

p  =  fluid  viscosity 

t  =  time 

UiyUj  =  x-component  of  fluid  velocity  component  tangential  to  lower  and 
upper  bearing  surface,  respectively. 

^1*^2  =  z-comPonent  of  fluid  velocity  component  tangential  to  lower  and 
upper  bearing  surface,  respectively. 

The  following  assumptions  are  inherent  in  the  development  of  Reynolds'  equation: 
o  The  lubricant  is  Newtonian  —  that  is,  the  shear  stress  in  the  lubricant 
is  proportional  to  the  rate  of  shear, 
o  The  flow  is  laminar. 

o  The  fluid  film  thickness  is  negligibly  small  in  comparison  with  the 
overall  scale  of  the  bearing, 
o  Fluid  inertia  is  negligible. 

o  The  lubricant  adheres  to  the  surfaces  (no  slip), 
o  No  external  forces  act  on  the  lubricant, 
o  The  lubricant  viscosity  is  uniform  across  the  film. 

A  thrust  bearing  usually  has  sector  pads.  It  is  more  natural  to  study  the  sector 
pad  in  cylindrical  polar  coordinates.  The  corresponding  Reynolds'  equation  is 


d  /  h3p  dp  \  d  /  rh3 p  dp  \ 

rdd  \  u  rdO  /  rdr  \  u  dr  / 


d  d 

=  6  2  -  (ph)  +  —  Iph  (oj.+coJI 
dt  dQ  12 

where:  r  =  the  radial  coordinate  in  polar  notation 

9  =  the  circumferential  angular  coordinate  in  polar  notation 

G>j,a>2  =  the  angular  rotational  speed  of  the  lower  and  upper  bearing 
surfaces 


In  allowing  (cOj.g^)  to  appear  as  independent  quantities,  a  duo-rotation  system  can 
be  considered.  Often  the  rectangular  form  approximation  is  used,  and  in  the  case 
of  an  incompressible  lubricant,  one  has 

<3  /  h  dp  \  d  /  h  dp  \ 
dCC  '  /i  dx  ■  dz  '  dz  / 


=  6 


dh  dh  d 

2—  +OJ+U —  +  h  —  iU+UA 
dt  1  2  dx  dr  1  - 


The  terms  on  the  right  hand  side  of  the  equation  are  referred  to  as  the  squeeze 
term,  the  wedge  term,  and  the  stretch  term,  respectively.  The  stretch  term  is 
seldom  of  any  importance.  This  form  of  equation  can  be  used  to  study  the  two- 
dimensional  tilt  pad  bearing  illustrated  in  Figure  2.  In  this  particular 
instance,  the  surface  of  the  tilt  pad  does  not  have  a  sliding  motion  (U^  =  0). 


FILM  VELOCITY 
DISTRIBUTION 


Figure  2  -  Two-Dimensional  Tilt  Pad  Bearing;  Top  View  Shows  Width, 

b,  and  Length,  L,  of  Bearing 


The  pressure  distribution  in  the  lubricant  fiLm  can  be  determined  by 
integrating  in  the  x  and  z  direction,  yielding  the  load  capacity,  FN,  as 

fn  =  f  f  dx  dz 

JO  h 

The  flow  per  unit  width  in  the  x  direction  is 


and  in  the  z  direction 


h  h  dp 

q  =  (U.  +  UJ  -  - - 

1  1  2  2  \2p  dx 


h*_  dp 
q*  ~  ~  12  dz 


The  friction  force,  Fp,  in  the  lubricating  film  is  the  sum  of  the  velocity-induced 
and  pressure-induced  shear  stresses  as 

r  L  rb  p{V2-Vy)  !  rc  rb 


Ff=± 


1  O 
-1 

2  o  o  dx 


J  0  '  0  h  2  Jo  Jo  dl 

where  the  plus  and  minus  signs  are  used  for  the  force  acting  on  the  lower  and 
upper  surfaces,  respectively.  The  film  thickness  must  be  describable  in  terms  of 
x  and  z  for  use  in  Reynolds'  equation. 

Other  simplifications  are  usually  possible  in  applying  Reynolds'  equation  to 
hydrodynamic  bearings.  For  steady  state  conditions  the  squeeze  term,  dh/dt, 
vanishes.  If  the  bearing  and  mating  surfaces  are  made  of  rigid  material,  the 
stretch  term  also  vanishes.  In  most  applications  the  bearing  speed,  Ui,  is  zero. 
Reynolds'  equation  now  becomes,  writing  U  =  U2 


±( 

dp  ' 

)+  Li 

(  h 3 

dp  \ 

dx  ' 

1  p 

dx  ) 

]  dzl 

P 

dz  / 

If  the  viscosity  can  be  considered  constant  over  the  entire  domain,  further 
simplification  can  be  made: 


±(*»*)+  =  6„U^ 

dx  V  dx  /  dz  [  dz  1  dx 


This  equation  is  most  commonly  used  to  determine  the  pressure  distribution,  from 
which  one  can  determine  the  bearing  load  capacity  and  the  friction  for  operation 
that  is  free  of  asperity  interactions  (hydrodynamic  conditions). 


Reynolds'  equation  has  also  been  applied  to  the  region  of  mixed  lubrication, 
where  load  is  supported  by  both  hydrodynamic  pressure  and  asperity  contacts. 
Christensen^  proposed  stochastic  models  to  permit  consideration  of  surface 
roughness.  For  the  one-dimensional,  longitudinal  roughness  (Figure  3),  he  recast 
the  Reynolds'  equation  in  terms  of  mean  or  expected  values  of  pressure  and  film 
thickness  parameters^: 

eM-eU  +  -  [— l—  -fifpH  =6p(/-Eff/) 


'ax  IT  lobe  Til  dz  l  £{1///3}  dz  ~  r  J  I  6flU  dx  Er\ 
where:  H  =  the  fluid  film  thickness,  including  roughness  deviations 
E{  }  =  the  expectation  or  stochastic  average  of  {  } 

The  corresponding  film  thickness  profile  was  treated  as 

H  =  h(x,z)  +  hs(z) 

Where  h  is  the  smooth  part  of  the  film  thickness  profile  and  hs(z)  represents  the 
roughness  deviations,  which  are  assumed  to  be  on  one  of  the  two  surfaces. 
Independence  of  hs  on  x  indicates  the  longitudinal  topography.  Stochastic 
description  of  the  roughness  deviations  is  represented  by  the  probability 
distribution  function,  which  is  assumed  to  have  the  polynomial  form  of 


P(/i  )  = 

S 


-  c<  h  <  c 

s 


where:  P(hs)  =  the  probability  distribution  function  of  hs 

c  =  half  the  total  range  of  the  random  film  thickness  variable,  S 


E  {}is  the  ensemble  average  with  respect  to  the  assumed  probability  distribution 
function.  As  a  model  problem,  Christensen  analyzed  in  detail  the  infinitely  long, 
tapered  slider. 

Christensen's  approach  was  applied  to  a  crowned,  tilt-pad  bearing  (shown  in 
Figure  A)  under  a  David  Taylor  Naval  Ship  R&D  Center  contract.  Crowning  was 
introduced  to  simulate  initial  pad  crown  and  distortions  during  operation  due  to 
pad  loading  and  thermal  gradients.  The  tilt-pad  feature  required  the  pad 
inclination  to  be  determined  by  the  condition  of  moment  equilibrium  about  the 
pivot  point.  This  condition  included  the  contribution  of  asperity  contact 
friction,  which  would  become  prominent  in  mixed  lubrication.  The  technical  report 
of  this  analysis  is  included  as  Appendix  A.  Although  the  predicted  frictional 
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performance  of  the  tilt-pad  thrust  bearing  in  the  hydrodynamic  and  mixed 
lubrication  region  resemble  experimental  results,  several  theoretical  issues 
remain.  These  can  be  grouped  as  follows: 
o  Film  striation  or  cavitation 
o  Thermo-elastohydrodynamic  effects 
o  Computation  techniques 
o  Surface  roughness  effects 

In  the  rest  of  this  report  we  have  addressed  theoretical  approaches  to  various 
topics  under  these  headings.  These  approaches  may  serve  as  guidelines  for  future 
in-depth  studies. 

THEORETICAL  TOPICS  RELATED  TO  THE  STUDY  OF  MIXED  LUBRICATION 
OF  SLIDING  SURFACE  BEARINGS 

1.  FILM  STRIATION 

Film  striation  often  occurs  in  sliding  bearings.  Part  of  the  bearing  gap  is 
filled  with  gaseous  and/or  vaporous  constituents  at  close  to  ambient  pressure. 

The  boundaries  of  the  striation  region  must  be  accurately  determined  for 
computation  of  the  film  pressure,  the  flow  rate,  the  viscous  shear  stress,  and 
viscous  heating.  The  following  derivations  incorporate  film  rupture  into  the 
computation  of  pressure  and  shear  in  a  fluid  film. 
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If  a  full  film  is  allowed  to  occupy  the  entire  shoe  length,  then  the  constant 
of  integration,  Ce,  is  to  satisfy  the  ambient  condition  imposed  at  the  exit  point 
of  the  fluid  film: 


p(b)  =  0 


(1-3) 


or 


_  E{H]dz  /\h 

Jo  £{/J3}/  'o 


dx 

E{H3} 


(1-4) 


as  indicated  in  Appendix  A.  However,  in  the  case  of  a  converging-diverging  film 
profile,  Equations  (1-2)  and  (1-4)  can  result  in  a  subambient  pressure  domain 
immediately  preceding  the  exit  edge.  The  latter  circumstance  violates  common 
experience.  It  is  generally  believed  that  a  subambient  pressure  condition  in  a 
fluid  film  cannot  be  sustained;  instead,  the  film  ruptures  or  cavitates  somewhere 
so  that  the  ambient  condition  will  be  reached  internally  with  a  vanishing  pressure 
gradient.  This  is  known  as  the  Reynolds-Swift-Stieber  condition  of  film  rupture 
or  cavitation  and  is  stated  as 


p(xc)  =  0 


e 


0  <  xc  <  b 


To  satisfy  Equation  (1-6),  Equation  (1-1)  yields 


CF  =  -E{H) 


(1-5) 


(1-6) 


(1-7) 


(1-8) 


Substitute  Equation  (1-8)  into  Equation  (1-2),  setting  x=xc;  one  obtains 


p(xc)  =  6/iU 


c  E{H]dx 
)  E{H3} 


-E{Hi 


I  c  _dr 
'o  Etf3} 


This  may  be  regarded  as  a  transcendental  equation,  the  root  of  which  yields  the 
location  of  the  beginning  of  film  rupture  or  the  cavitation  domain. 

The  necessary  steps  to  determine  the  root  of  Equation  (1-9)  can  be  readily 
incorporated  into  a  computer  program.  Equation  (1-9)  can  be  regarded  as  a 
residual  function  and  can  be  obtained  within  the  cycle  of  numerical  quadratures 
for  the  two  integrals.  Sign  reversal  of  the  residual  function  signifies  that  the 
beginning  of  film  rupture  had  just  been  overtaken,  and  a  simple  interpolation 
would  yield  its  location  and  the  correct  value  of  Cg  quite  accurately. 

In  common  bearing  configurations,  the  bearing  gap  continues  to  diverge  until 
the  exit  edge  of  the  bearing  pad  is  reached,  so  that  there  is  no  need  to  consider 
whether  the  full-film  condition  can  be  re-established.  Even  the  full  cylindrical 
journal  bearing  usually  has  a  feed-groove  or  a  feed-hole  located  in  the  unloaded 
part  of  the  bearing  surface,  so  that  it  may  be  regarded  as  a  circular  arc  of 
approximately  360°,  and  the  maximum  gap  point  can  be  assumed  to  be  the  end  of  the 
rupture  domain  or  the  full-film  reestablishment  point. 


Viscous  Shear  of  a  Striated  Film 


The  film  flux  of  a  bearing  with  ruptured  film  is  ^UE(H)|Xc.  In  the  domain 
xc<x^b,  this  film  flux  is  insufficient  to  fill  the  local  film  thickness  E(H)|X. 
According  to  the  striation  model  of  film  rupture,  the  local  striation  length 
fraction  is  E(H)|X_/E(H) .  Accordingly  the  local  viscous  shear  stress,  t,  is 


The  analysis  presented  in  Appendix  A  is  a  closed  integral  solution  to 
Reynolds'  equation.  End  leakage  is  neglected.  During  operation  lubricant  flow 
occurs  not  only  in  the  direction  of  motion  (leading  edge  to  trailing  edge)  but 
also  leaks  along  ends  of  the  pad.  Such  end  leakage  can  significantly  affect  the 
load  capacity  of  the  bearing  by  bleeding  off  pressures  built  up  in  the  film. 
Fuller^  presented  an  approximate  procedure  to  account  for  end  leakage  through  the 
use  of  a  correction  factor  for  the  load  and  the  minimum  film  thickness.  The 
correction  factor  was  based  on  the  experimental  data  for  a  tapered  wedge  with 
bearing  width-to-length  ratio  ranging  from  ^  to  5J.  Such  an  approximate  procedure 
for  correcting  end  leakage  can  be  used  in  conjunction  with  the  analysis  of 
Appendix  A.  A  later  section  of  this  report,  Finite  Difference  Computation  of  the 
Infinitely  Long  Slider  (in  the  section  on  computation  techniques  under  "High 
Accuracy  Numerical  Algorithm  for  Computing  Fluid  Film  Pressure"),  presents 
approaches  for  considering  the  use  of  finite  difference  techniques  to  prepare  a 
more  thorough  treatment  of  end  leakage. 

2.  THERMO-ELASTOHYDRODYNAMIC  EFFECTS 

General  Discussion 

Thin-film  lubrication  is  prevalent  in  many  critical  machine  elements  in  naval 
propulsion  and  shipboard  machine  systems.  Trouble-free  operation  of  such  machine 
elements  depends  on  the  control  of  subtle  geometrical  details  of  the  film,  details 
that  are  readily  influenced  by  elastic  and  thermoelastic  deformations  of  the 
surfaces.  Rational  design  procedures  for  such  machine  elements  have  not  been 
developed.  The  derivations  in  the  following  sections  provide  general  guidance  for 
the  development  of  such  analytical  engineering  tools. 

The  possibility  of  adjusting  surface  geometry  for  optimization  of  designs  was 
eloquently  demonstrated  by  Rayleigh's  classic  paper^,  and  elastohydrodynamic  and 
thermo-elastohydrodynamic  effects  in  thin-film  machine  elements  are  especially 
important  in  naval  applications. 


For  various  reasons,  lubrication  by  very  thin  fluid  films  in  naval  machine 
elements  is  likely  to  be  associated  with  pressure  distribution  that  can 
significantly  deflect  the  lubricated  surface.  At  the  same  time,  the  viscous  shear 
stress  in  the  fluid  film,  which  is  responsible  for  the  generation  of  load-carrying 
film  pressure,  necessarily  also  generates  heat  and  accompanying  temperature  rise 
and  temperature  gradients  in  the  machine  parts.  Important  thermal  phenomena 
include  reduction  of  fluid  viscosity,  differential  thermal  growths  (due  to 
different  coefficients  of  thermal  expansion  of  various  materials  in  the  machine 
structure),  and  thermal  distortions  (due  to  temperature  gradients  regardless  of 
the  temperature  level). Three  major  categories  of  machine  elements  deserve 
particular  recognition  in  relation  to  studies  of  thermo-elastohydrodynamic 
effects.  An  example  in  the  first  category  is  the  thrust  bearing  in  a  ship 
propulsion  system,  which  carries  the  load  caused  by  sea  pressure  acting  on  an 
unbalanced  area  and  by  vessel  acceleration.  The  study  of  tilting-shoe  thrust 
bearings  first  revealed  the  significance  of  thermo-elastohydrodynamic  effects. 
Thermally  induced  surface  crowning  is  probably  responsible  for  the  stabilization 
of  the  pitch-roll  motions  of  the  bearing  shoes.  An  example  in  the  second  category 
is  the  dynamic  seal,  which  not  only  keeps  water  out  of  the  ship  interior  but  also 
prevents  lubricants  and  process  fluids  from  contaminating  the  cabin  environment. 
Although  "contact"  operation  is  often  assumed  for  low-leakage  seals,  current 
belief  holds  that  low  wear  rate  is  critically  dependent  on  the  effective 
generation  of  film  pressure  to  maintain  a  stable,  partial  separation  of  the  mating 
surfaces.  Several  mechanisms  have  been  suggested  to  explain  the  likely  existence 
of  load-carrying  film  pressure  in  low-leakage  face  seals:  namely,  hydrodynamic 
wedge  action  due  to  runner  misalignment  or  thermoelastic  waviness  of  the  mating 
surfaces;  pressure-induced  gap  convergence;  and  evaporation  of  leakage  fluid.  An 
example  from  the  third  category  of  machine  elements,  one  that  has  fascinated 
researchers  studying  probable  thermo-elastohydrodynamic  effects,  is  the  sterntube 
bearing,  which  is  constructed  of  elastomeric  staves.  A  relatively  large 
coefficient  of  thermal  expansion  and  a  relatively  small  shear  modulus  make  the 
elastomeric  stave  bearing  a  prime  subject  for  study  of  both  the  beneficial  and  the 
harmful  effects  that  have  been  attributed  to  these  properties.  Table  1 
illustrates  the  essential  thermo-elastohydrodynamic  phenomena  associated  with 
these  three  machine  elements. 
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The  generic  machine  element  of  interest  consists  of  two  bodies,  which  are 
respectively  designated  as  (1)  and  (2),  and  a  thin  fluid  film  that  permits  low 
friction  relative  sliding  between  the  two  bodies.  The  internal  mutual  action 
between  the  two  bodies  is  transmitted,  for  the  most  part,  through  the  enveloping 
surfaces  of  the  thin  film  in  the  form  of  distributed  traction.  Friction  heating 
in  the  thin  film  is  due  to  viscous  dissipation  and  may  also  be  associated  with 
asperity  contact  sliding.  The  latter  can  precipitate  tribological  trauma  and  thus 
cause  catastrophic  failure.  Under  ordinary  circumstances,  friction  heating  is 
converted  into  sensitive  heat,  raising  the  temperature  of  the  two  bodies  and  the 
circulating  fluid. 

Conservation  laws  of  mass,  momentum,  and  energy  as  applied  to  each  of  the 

three  parts  of  the  system  and  as  applied  to  the  overall  system  are  the  foundation 

for  the  thermo-elastohydrodynamic  analysis  of  machine  elements.  The  essential 
elements  of  the  mathematical  statements  of  this  problem  are  derived  below. 

Conservation  Laws  for  the  Fluid  Film 

The  significance  of  thermo-elastohydrodynamic  analysis  for  machine  elements 
was  broadly  discussed  above  in  connection  with  naval  machine  systems.  The  generic 
machine  element  is  a  system  consisting  of  two  bodies  separated  by  a  thin  film  of  a 
viscous  fluid.  The  functional  behavior  of  the  machine  element  depends  on  the 
dynamic  and  thermal  interactions  between  the  thin  fluid  film  and  the  two 

surrounding  bodies.  The  present  section  summarizes  the  applicable  equations 

specifically  pertaining  to  the  thin  fluid  film. 

For  an  incompressible  Newtonian  lubricant,  the  mathematical  theory  of 
lubrication  is  usually  presented  in  terms  of  the  Reynolds'  equation,  which  is 
summarized  in  the  contemporary  form  as  follows: 


where:  V  =  the  film  flux  vector 


Jq  =  the  integral  of  the  viscosity  weighted  Couette  velocity  profile 
function 

Jp  =  the  integral  of  the  viscosity  weighted  Poiseuille  velocity  prof 
function 

V  =  che  two-dimensional  surface  gradient  operator 

(J>  =  leakage  flux  through  bearing  surfaces 

In  the  "surface"  vector  form  presented  here,  these  expressions  are  valid  for 
flat  and  moderately  curved  surfaces.  Figure  5  illustrates  a  generic 
representation  of  such  a  pair  of  surfaces.  The  film  geometry  has  two  length 
scales;  R  is  the  size  scale  of  the  bearing  while  C  is  the  separation  scale 
representing  the  nominal  bearing  clearance. 


Figure  5  -  Thin  Film  Between  Nearly  Parallel  Surfaces 


A  suitably  chosen  reference  surface  is  used  to  define  the  geometry  of  lower 
and  upper  bearing  surfaces  in  terms  of  their  respective  normal  distances  (h^,  l^). 
A  unit  normal,  n*,  can  be  defined  at  every  point  of  the  reference  surface;  it  can 
be  regarded  to  be  also  perpendicular  to  either  bearing  surface  and  is  usually 
"outwardly"  directed  by  convention.  The  three-dimensional  fluid  velocity  V  is 
described  by  its  normal  component  V*n  and  its  tangential  component  U;  the  latter 
can  further  be  decomposed  into  appropriate  surface  components  according  to  the 
representation  of  the  reference  surface.  ^ is  the  surface  gradient  operator, 
which  is  constrained  to  exclude  the  normal  component  of  the  conventional  three- 
dimensional  gradient  operator. 


The  radius  of  curvature  of  the  reference  surface  is,  by  definition,  of  the 
same  order  as  the  size  scale.  The  separation  to  size  ratio  C/L,  is  typically  10-^ 
to  10  .  A  necessary,  but  sometimes  overlooked,  assumption  is 


W(hvh.,) 


(2-2) 


Where  0{  }means  "order  of  {  }."  Any  discontinuities  in  (hj,  J^)  will  require 
special  attention.  M1  is  the  tangential  (to  the  reference  surface)  flux  vector  and 
is  obtained  by  integrating  the  tangential  velocity  component  with  respect  to  the 
separation  between  the  two  surfaces: 


(2-3) 


Subscripts  (1,  2)  as  applied  to  U  and  h  refer  to  the  lower  and  upper  surfaces, 
respectively.  The  film  pressure  p  is  regarded  as  invariant  with  respect  to  the 
film  thickness  coordinate  y  in  the  lubrication  theory,  h  =  l^-h^  is  the  local 
separation  between  the  surfaces.  £$  accounts  for  porosity  of  the  bearing  surfaces 
and  is  the  sura  of  the  (wall)  penetration  flow  rate  per  unit  area  of  the  two 


Simplification  of  the  Navier-Stokes '  equations  to  the  lubrication  theory,  as 
represented  by  Equations  (2-la)  and  (2-lb),  requires  a  number  of  prior 
assumptions.  Let  V  be  the  scale  of  fluid  velocity.  The  necessary  assumptions  to 
ensure  validity  of  the  lubrication  theory  include  the  following: 

(1)  Components  of  the  dyad  V  V  are  dominated  by  the  term  Ctf'V )  U  which  is  of  0 

{V/C}.  Consequently,  the  stress  tensor  enters  into  the  momentum  conservation 
law  through  the  term 

V-I  =  -Vp  ♦  Cfrv)  IpCrf-VHf) 

where:  V  =  the  three-dimensional  gradient  operator 

V*i  =  the  dot  product  of  the  three-dimensional  gradient  operator  and  the 
viscous  shear  stress 


=  the  unit  vector  perpendicular  to  the  nominal  bearing  surface 
=  the  fluid  velocity  vector 

p0  is  the  nominal  viscosity  coefficient;  the  actual  viscosity 
coefficient  may  vary  throughout  the  lubricant  film. 


(2)  The  film  Reynolds'  number  (pVC/p0)  (C/L)  is  of  negligible  magnitude. 

(3)  Mass  conservation  is  satisfied  with  V‘7  and  (n*\^)  (rrV^),  cancelling  out  each 
other  so  that 


-»  -»-»  CV  1 

\U\=0{V\-  |n-V|  =  0  — 

(4)  The  flow  field  is  laminar.  Consistent  with  these  assumptions,  the  fluid 

pressure  can  be  shown  to  have  negligible  variation  in  the  normal  direction. 
The  momentum  conservation  law  is  reduced  to 

6  (  dU 
dy  ■  dy 

(2-4 
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which  can  be  integrated  with  respect  to  y  twice;  upon  imposing  nonslip 
conditions  at  the  surfaces,  one  obtains 

7=  Uj  +  df2  -  1^)  Ic  -  (Vp)lP 


where:  Ic  =  /g(y)//Q(h2) 

Ip  =  /j(h2)  IC  -  I±( y) 


(2-5) 


'  »w  =  \ 


y  ykdy' 


(2-6a,b,c) 


1^  is  the  kth  moment  integral  of  inverse  viscosity.  Iq  and  Ip  are,  respectively, 
velocity  profile  functions  associated  with  shear  and  pressure  effects. 

The  lubricant  viscosity,  at  common  operating  environments  of  fluid  film  bearings, 
should  be  regarded  as  temperature  dependent.  If,  however,  one  chooses  to  accept 
the  isoviscous  approximation,  p  =  pQ  =  const.,  then 

(y-V 

^ c\soinscous  (h  —h) 


<IJ 

r  isovisants 


lh2-y)(y-  h{) 


(2-7a,b) 


These  functions  are  seen  to  describe  the  Couette  and  Poiseuille  velocity  profiles, 
respectively.  The  coefficients  (Jq,  Jp)  in  Equation  (2-la)  are  obtained  by 
substituting  Equations  (2-5)  and  (2-6a,  b)  into  Equation  (2-3)  and  carrying  out 
suitable  integrations  formally: 


(2-8a,b) 


JC  |  lCd*  JP 
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And  for  the  isoviscous  approximation, 


(JJ  =  ~  (A  -  h.)  = 

C  ISOVISCOUS  2  2  * 


1 

-A 

2 


(JJ 

r  isoviscous 


(*2  -  V 


12^f 


12/ir 


(2-9a,b) 


Accurate  determination  of  the  film  pressure  requires  a  suitable  allowance  for 
the  temperature  dependence  of  viscosity,  so  that  Equation  (2-8a,b)  instead  of  (2- 
9a, b)  can  be  used.  Therefore,  it  is  of  interest  to  examine  the  energy 
conservation  law  as  applied  to  the  lubricant  film  flow  process.  The  complete 
governing  equation  of  energy  in  a  liquid  is 


DE 


DT 


p  -  +  Bp -  =  Q  +  V-  ( kVT) 

H  Dt  Dt 


(2-10) 


where  E0  is  the  specific  internal  energy,  (5  is  the  bulk,  coefficient  of  thermal 
expansion,  Q  is  the  internal  dissipation,  k  is  the  thermal  conductivity,  and  T  is 
the  fluid  temperature.  (D/Dt)  is  the  substantive  time  derivative  which  includes 
the  convective  derivative  to  account  for  fluid  motion.  For  common  conditions  in 
fluid  film  bearings,  the  second  term  on  the  left  hand  side  of  Equation  (2-10)  can 
be  neglected.  Also,  consistent  with  the  approximations  already  invoked  in  the 
lubrication  theory,  the  terms  on  the  right  hand  side  may  be  simplified  as 


(2-lla,b) 


If  asperity  contacts  occur  because  of  surface  roughness,  the  dissipation  term, 
Equation  (2-lla),  should  be  increased  to  reflect  sliding  friction  effects. 


Also,  writing  DEg  =  CVDT,  Cy  being  the  specific  heat  for  a  constant  volume 
process,  then  the  applicable  energy  equation  which  defines  the  temperature  field 


in  the  fluid  film  is 


pC  —  =p(  —  )•(  —  )+  —  (*— ) 
v  Dt  V  dy  )  \  dy  /  dy  '  dy  / 


For  the  steady-state  problem, 


pC  (V-V)T 


(2-12 


.„(£).(£)♦*(.* I) 

\  dy  /  \  dy  /  dy\  dy  J 


(2-13 


This  equation  is  of  the  parabolic  type  and  its  rigorous  field  solution  is  to 
satisfy  an  inlet  temperature  profile  and  suitable  continuity  conditions  of  heat 
flux  at  the  surfaces.  Equation  (2-13)  can  also  be  expressed  in  the  thin  film 
integral  form: 

__  f  *2 

pC  W-  ( UT)dy 


(2-14 


Equation  (2-14)  can  be  utilized  to  seek  an  approximate  temperature  field  without 
dealing  with  the  full  complexity  of  Equation  (2-13).  Both  the  convective  term  on 
the  left  hand  side  and  the  dissipation  term,  the  first  term  on  right  hand  side, 
need  to  be  known  as  integrals  across  the  surface  separation;  they  can  tolerate 
crude  approximations  without  introducing  very  serious  errors.  The  last  term  is  to 
be  matched  with  conduction  heat  flux  entering  the  bodies  which  surround  the  fluid 


The  cylindrical  polar  coordinate  system  can  be  used  for  the  three- 
dimensional,  geometrical  representation  of  either  the  cylindrical  journal  bearing 
or  the  thrust  bearing.  The  complete  three-dimensional  gradient  operator  is 

-*d  -»  l  a  -»  a 

V  =  i  —  +  i  -  —  -hi  — 
r  dr  r  dd  2  dz 

The  components  of  a  vector  A  are  designated  by  the  subscripts  (r,0,z) 
respectively.  The  divergence  operation  takes  on  the  form 

-»  1  d  1  dAe  dAz 

V-A  = - (rA  )  + - +  - 

r  dr  r  r  dQ  dz 

The  sliding  velocity  is  due  to  rotation: 

U  =  igCjr 

If  both  surfaces  slide,  then  subscripts  (1,2)  would  be  associated  with  the 
rotational  speed,  o.  In  the  case  of  the  counter-rotating  shafting,  (o>j,  would 

have  opposite  signs. 

For  the  cylindrical  journal  bearing,  the  nominal  reference  surface  is  a 
circular  cylinder.  The  normal  direction  is  radial.  Equations  (2-la,b)  thus 


become 


-*-*[[  '7  dp  1  -*(  dp  \ 

qj=l0  [  |"lA  +  (W2  -  M\)JC  r"  7  5)1  +  U  "  to  ) 

Id  •>  dp  |  d  (  dp  \  dh 

O  =  —  —  o  h.  +  ( o „  -  co  )  Jr  r  -  J  }  -  —  (  J  —  )  +  —  +  E<t> 

r2  dd  1  2  1  c  p  del  d z [  p  dz  /  dt 


(2-15a,b) 


The  radial  coordinate  r  can,  within  the  accuracy  of  the  lubrication  theory,  be 
regarded  as  a  constant  and  be  equated  to  the  nominal  journal  radius,  R. 

For  the  thrust  bearing,  the  nominal  reference  surface  is  a  radial  plane.  The 
normal  direction  is  axial.  Reynolds'  equation  takes  on  the  form 

-*-*[[  I  7  dp  -*(  dp  \ 

v=iel  a\h  +  (w2  “  wi)7|r~  7  #  +'r(-JpJr  j 


23 


(2-16a,b) 


O  =  -  ^  (  L/* 

r2  <30  ll  1 


+  (co2-a>l)Jc 


2  ,  dP 


1  <3  /  dp\  dh 
-  -  -  r  J  -  +  —  +  Z4> 

r  dr  \  p  dr  /  St 


In  contrast  to  Equation  (2-15),  in  which  r=R  is  a  constant,  r  in  Equation  (2-16) 
is  the  radial  coordinate. 


To  restate  the  energy  equation  in  a  more  usable  form,  each  term  in  Equation 
(2-14)  will  be  separately  examined.  For  either  bearing  type,  the  integrated 
dissipation  can  be  rewritten  as 


2  0 

(o>.2  -  w)  r" 


Q  dy  = 


W 


+ 


Wp 


w 


(2-17) 


7q  and  I ^  are  moment  integrals,  respectively  of  the  Oth  and  2nd  orders,  of  the 
inverse  of  the  viscosity  coefficient  as  given  previously  by  Equation  (2-6c). 

The  surface  conduction  terms  do  not  need  further  attention.  The  convective 
heat  transfer  terms  differ  slightly  in  form  for  the  two  bearing  types,  mainly  in 
the  divergence  operator: 

f*2  -> 

pCvV-  (UT)  dy 


=  pC 


<30 


o).  0  +  {(o  —  (o  )  K 

1  2  1  r 


-W-mpWp) 


(2-18) 


where , 


(2-19a,b,c) 


For  the  cylindrical  journal  bearing, 

75(^5) 


while  for  the  thrust  bearing, 

W-  (KPWP)  =  ~7  ~  (kp  ~  )  +  ~  -  (xpr  ~  ) 

p  p-  86  '  p  dd  !  r  dr  V  p  dr  ) 


(2-20) 


(2-21) 


The  fluid  film  is  not  a  closed  system  regarding  the  conservation  laws.  The 
following  are  sources  of  interactions  between  the  fluid  film  and  the  two 
surrounding  bodies: 

(1)  Global  equilibrium  of  the  machine  element  requires  that  the  internal  load  and 
moment  transmitted  through  the  fluid  film  be  related  to  the  kinematics  of  the 
two  surrounding  bodies  and  the  externally  imposed  constraints. 

(2)  The  surface  separation,  h  =  l^-hp  which  is  the  dominant  factor  in  the 
generation  of  the  film  pressure,  is  significantly  affected  by  thermoelastic 
deformations  of  the  two  bodies.  Thermoelastic  deformations,  in  turn,  result 
from  temperature  fields  caused  by  viscous  heating  which  is  moderated  by 
convective  and  conductive  cooling  effects  as  described  by  Equations  (2- 
14,17,18).  A  major  computational  difficulty  is  associated  with  the  accurate 
determination  of  the  temperature  field  with  its  simultaneous  conduction- 
convection  effects. 
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(3)  The  surface  separation  is  also  affected  by  elastic  deformations  caused  by 

film  pressure.  This  aspect  is  quite  important  in  sterntube  bearings  because 
of  the  compliant  elastomeric  staves.  It  is  also  important  for  high-pressure 
face  seals.  Load-induced  cross-sectional  twisting  of  the  seal  ring  can  be 
considerable  in  comparison  with  the  thickness  of  the  film  formed  by  leakage 
flow.  Elastic  deformation  from  film  pressure  is  less  important  for  tilting- 
shoe  thrust  and  tilting-shoe  journal  bearings,  where  pitch-roll  equilibrium 
becomes  an  inherent  feature  of  the  overall  analysis. 

(A)  Temperature  distribution  across  the  film  thickness,  in  addition  to  the 
surface  temperature  fields,  influence  the  film  pressure  because  of  the 
temperature  dependence  of  lubricant  viscosity.  Equations  (2-6a,b,c)  pertain 
to  this  issue. 

(5)  Surface  permeability  due  to  material  porosity  can  also  affect  film  pressure. 

It  can  be  important  in  the  friction  behavior  of  the  retainer-roller  contact 
of  rolling  element  bearings. 

Thermal  Effects 

Thermal  coupling  between  the  fluid  film  and  surrounding  bodies  stems  from  two 
independent  physical  causes:  (1)  temperature  dependence  of  lubricant;  and  (2) 
thermal  distortion  of  bearing  surfaces. 

The  Walther-ASTM  formula  may  be  used  to  determine  the  kinematic  viscosity,  v, 
from  two  data  points: 

log  log  (v  +  0.6)  =  log  T  +  m2 

(2-22) 

where  v  and  T  are,  respectively,  in  units  of  centistokes  and  degrees  Rankine  and 
m^  and  m2  are  coefficients  which  depend  upon  the  oil.  If  the  temperature 
distribution  within  the  fluid  film  is  known,  then  the  weighted  integrals  of  the 
reciprocal  of  viscosity  coefficient  can  be  accordingly  calculated: 


*  yV 

h{  niT) 


where  y'  is  a  dummy  variable  of  integration  for  y 


ic(y)  =  /0CyV/0(V 


(2-23) 


lpiy)  =  l]{h2)lc{y)-Il(y) 


J 

C 


lcdy 


(2-24a,b) 


J  = 

P 


(2-25a,b) 


Subsequently,  the  film  pressure  can  be  found  by  solving  the  generalized 
Reynolds'  equation  for  a  fluid  contained  between  nonpermeable,  smooth  walls: 

tnj  d/l 

v-v  +  —  =  o 

dt 


V  =  U{h  +  ((/2  -  UX)JC  -  Ctyp)jp 

(2-26a,b,c ) 

— ¥  — > 

Uj  and  U2  are  sliding  velocity  vectors  of  the  lower  and  upper  surfaces, 
respectively.  The  film  velocity  profile  is 

u=ux  +  (u2-ux)\c-(Wp)\p 

(2-27) 
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The  temperature  field  itself  is  to  be  found  by  solving  the  thin-film  energy 
equation: 


pC  ( V-V)T  = 

r  V 


\  dy  J  \  dy  /  dy\dy/ 


(2-28) 


V  is  the  three-dimensional  velocity  vector. 
Equation  (2-27),  yielding 


3U 

dy 


f'W 


-  (Wp)( 


dU/dy  can  be  found  by  differentiating 


w 

w 


(2-29) 


Equation  (2-28)  is  truly  three-dimensional;  it  is  of  the  parabolic  type,  being 
first  order  in  the  tangential  direction  and  second  order  in  the  normal  direction. 
Needless  to  say,  its  complete  solution  is  a  formidable  task.  Thus,  it  is  often 
desirable  to  accept  the  compromise  of  tackling  the  integral  form  of  the  energy 
equation : 


pCW- 


( VT)dy  = 


W 


Wp 


w- 


w 


w 


dT 

+  K  - 


(2-30) 


Analysis  of  the  Temperature  Field 

The  fluid-film  is  coupled  to  the  two  surrounding  bodies  through  the 
equilibrium  conditions  (both  static  and  dynamic),  elasticity  effects,  and  thermal 
coupling.  The  thermal  coupling  problem  will  be  treated  here. 

Consider  the  schematic  illustration  given  in  Figure  6.  Bodies  I  and  II  are 
separated  by  the  thin  fluid  film.  Flow  through  the  thin  film  is  maintained  by  a 
supply  through  a  manifold  arrangement  in  Body  I.  The  exposed  surfaces  of  Body  I, 
marked  "S,"  are  submerged  in  a  fluid  bath  that  is  maintained  at  a  constant  ambient 
temperature.  The  thermal  environment  is  similar  for  Body  II.  Subsequent 
discussion  will  refer  specifically  to  Body  I  only. 
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BODY  II 


'THIN  FLUID  FILM 


Figure  6  -  Schematic  Thermal  System 


The  steady-state  temperature  field  in  Body  I  would  obey  the  classical 
potential  theory: 

?  =  -k,VT 
V-~q  =  .  V -  f  /r/  VT)  =  0 


(2-31a,b) 


where:  q  =  the  conductor  heat  flux 


/Cj  =  the  thermal  conductivity  of  Body  I 

If  the  body  is  moving,  for  example  the  runner  of  a  thrust  bearing,  the 
divergence  law  should  be  replaced  by  the  diffusion  law,  which  allows  for  the 
contribution  of  convection.  Nevertheless,  because  of  the  large  thermal  capacity 
of  the  moving  material,  an  isothermal  condition  is  approximately  maintained  along 
its  interface  with  the  fluid  film.  Thus  Equation  (2-31b)  in  a  two-dimensional 
form,  would  remain  applicable  even  to  the  runner.  The  coefficient  of  thermal 
conductivity  should  be  spatially  uniform;  consequently,  the  divergence  law  is 
equivalent  to  the  Laplace's  equation  V^T  =  0. 


Three  types  of  thermal  boundary  conditions  are  to  be  imposed: 
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(1)  At  the  surface  of  the  thin  fluid  film  there  is  a  yet  unknown 
distribution  of  normal  heat  flux. 

(2)  At  each  of  the  exposed  surfaces,  a  suitable  heat  transfer  coefficient 
can  be  assigned  to  relate  the  local  surface  temperature  rise  above  the 
ambient  to  the  heat  flux: 


-K~n-(VT)  =  hT(Ts-TA) 


(2-32) 


where:  h^.  =  the  heat  transfer  coefficient  at  the  exposed  surface 
Tg  =  the  temperature  at  the  exposed  surfaces  S 
Ta  =  the  temperature  of  the  ambient  fluid 

(3)  At  the  walls  of  the  supply  manifold,  an  adiabatic  condition  could  be 
assumed.  On  an  actual  bearing,  however,  a  feed  groove,  instead  of  a 
closed  manifold,  is  more  typical.  Because  circulation  between  the  feed 
flow  and  the  ambient  fluid  is  not  inhibited,  part  of  the  overall  cooling 
would  take  place  at  the  walls  of  the  feed  groove.  The  fluid  in  the  feed 
groove  should  be  somewhat  hotter  than  the  ambient  fluid.  Thus,  the 
adiabatic  supply  model  would  overestimate  the  temperature  rise  of  Body 
I.  The  other  extreme  model  would  assume  that  an  ambient  condition 
prevails  in  the  feed  groove  and  would  select  a  suitable  heat  transfer 
coefficient  for  the  feed  groove.  This  model  thus  would  underestimate 
the  temperature  rise  of  Body  I.  These  two  extreme  models  may  be 
regarded  as  upper  and  lower  bounds  of  the  thermal  situation. 

Standard  methods  are  available  for  solving  the  conduction  problem.  The 
boundary  condition  at  the  exposed  surfaces,  as  given  by  Equation  (2-32),  can 
accommodate  a  full  range  of  possibilities,  which  include  the  extremes  of  a  fully 
chilled  condition 


Ts  -  TAatS  if  hjJKj  -+<*> 
and  also  the  insulated  condition 

~n-(VT)-OatS  ifh^Kf-*0 


=  -  k  n  •  (VT)j-  at  the  interface. 

Here,  the  subscript  "i"  refers  to  body  I  and  the  subscript  "f"  refers  to  the  fluid 
film.  The  temperature  also  should  be  continuous  at  the  interface;  or 

(T)r  =  (T)^  at  the  interface. 

If  Body  I  is  moving  and  is  also  well  chilled,  then  (T)I  would  be  a  known 
constant  and  the  conduction  problem  of  Body  I  merely  serves  to  impose  a 
temperature  boundary  condition  at  the  interface  with  the  fluid  film.  It  is  no 
longer  necessary  to  treat  the  conduction  problem.  The  heat  flux  across  the 
interface  can  be  calculated  only  from  the  temperature  field  in  the  fluid  film. 

To  illustrate  the  relationship  between  the  temperature  and  the  heat  flux  at 
the  interface  with  the  fluid,  consider  a  model  two-dimensional  problem  for  the 
domain  (0  ^  X  SI  b,  0  S  yj  s  w).  Discretization  is  accomplished  by  dividing  the 
sides  a  and  b  respectively  into  M  and  N  uniform  intervals.  The  mesh  index  i  has 
the  range  of  (1  to  M+l)  in  the  x-direction  while  j  has  the  range  of  (1  to  N+l)  in 
the  yj-direction.  Let  j=l  be  the  interface  with  the  fluid  film,  while  the 
convectiveiy  cooled  boundary  condition  is  imposed  for  the  sides  (i=l,  i=M+l,  and 
j=N+l).  The  entire  temperature  field  is  represented  by  the  column  vectors 

(T}j 

for  ( j=l ,2 ,3, . . .N+l ) .  Each  vector  has  a  rank  of  M+l  and  contains  the  elements 
Tj^j  for  (i=l,2,3, . . .M+l) .  The  conduction  heat  flux  at  the  interface  directed 
into  the  fluid  film  is  the  vector 


which  is  also  of  rank  M+l.  The  conduction  problem  of  Body  I  can  thus  be  described 
in  terms  of  the  following  system  of  matrix  equations: 

[ A(°> ]t  {T},+[A<  +  >]t{T}2  =  .  {q}. 


(2-34) 


[A(  )]}  {t}j  j  +  [A<o)]j  (t}j  +  [a(+)]j  {t}j+j  =  {0} 

for  j=2,3,...N 

/■a,j;n+,  mn + [a'°’]n+,  =  {0} 


(2-35) 


This  system  can  be  solved  by  the  method  of  forward-elimination-and-backward- 
substitution,  resulting  in 


ffl,  =[B[]{qJ 

In  turn,  this  last  expression  can  be  inverted  to  yield 

M=[bi]i{t}1 


(2-36) 


(2-37) 


In  this  manner,  the  conduction  problem  is  reduced  to  a  pair  of  reciprocal 
influence  matrices  [Bj]  and  which  relate  the  interface  temperature  (T}j 

linearly  with  the  heat  flux  {q}j  and  vice  versa  according  to  Equations  (2-36)  and 
(2-37). 

Next,  the  energy  equation  of  the  fluid  film  will  be  treated  for  the  domain 

(0  S  x  ^  a,  0  5  yf  ^  h).  Equation  (2-28)  can  be  written  as 
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(2-38) 

The  subscript  "f"  has  been  dropped  from  yf  to  simplify  notation.  The  velocity 
components  (u,v)  are  regarded  to  be  given.  The  interface  with  Body  I,  temperature 
and  heat  flux  are  continuous.  Thus 


T(x,h)  =  T {x) 

Tj(x)  in  the  discretized  representation  is  simply  {T}^.  And 


(2-39) 
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Note  that  the  sign  in  Equation  (2-40)  is  consistent  with  the  definition  of  qj, 
which  is  directed  into  the  fluid  fiLm.  In  the  above  pair  of  equations,  neither  Tj 
nor  qj  is  fully  specified.  However,  they  are  related  to  each  other  according  to 
Equation  (2-36)  or  (2-37).  To  complete  the  physical  description  of  the  model 
problem  let  the  inlet  fluid-film  temperature  be  a  known  constant 

T(x=0.y)  =  Tin[e( 

(2-41) 

and  the  surface  temperature  of  Body  II  be  also  a  given  constant  to  represent  a 
runner  face: 


T  (x,y~0) 


T 

runner 


(2-42) 


Although  the  above  statements  are  sufficient  to  allow  computation  of  the 
temperature  field  in  the  fluid  film,  one  may  observe  that  Equation  (2-38),  which 
contains  the  convective  terms  on  the  left-hand  side  and  conduction-dissipation 
terms  on  the  right-hand  side,  is  a  troublesome  partial  differential  equation  from 
the  standpoint  of  numerical  analysis.  A  special  scheme,  which  consists  of 
iterative  successive  approximations,  is  proposed  to  simplify  the  method  of 
computation.  Essential  steps  of  the  proposed  iterative  scheme  are  outlined  below 
for  the  model  problem. 

(1)  Postulate  the  transverse  temperature  profile  to  be  of  the  form 


tv-  vj  =  t  -h  (T r  •  T  )  F  (x,y/ h) 

iix.yi  —  i  runner  l  runner 


(2-43) 


to  satisfy  Equations  (3-39)  and  (2-42),  profile  function  F(x,y/h)  should  be 


unity  at  y=h  and  vanish  at  y=0.  To  start  the  iterative  computation,  let 
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F(x,y/h)  —  y/h 


(2-44) 


for  the  first  time. 

(2)  Consider  the  integral  form  of  Equation  (2-38)  with  the  profile  function: 

'c.sCH;)r.K-r.(( f)  -(f)  | 


=  -pC  j  I  \u(l-F)T  dy-KT  (—)  -(  —  ) 

u  dx  J0  l  runner  runner  IV  dy  )ysh  \dy/=0l 


V2  (dp\  f  V 


(2-45) 


Note  that  this  is  an  ordinary  differential  equation  for  Tj  (x). 

(3)  Having  obtained  Tj  (x),  use  Equation  (2-37)  to  compute  qj,  which  then  yields 
icd T/8y  at  y=h. 

(4)  Integrate  Equation  (2-38)  across,  beginning  from  y=h,  with 


(-) 

'  %  >y  =  h  «■ 


(2-46) 


as  a  boundary  condition.  For  the  terms  in  the  left  hand  side  of  Equation  (2- 
38),  the  profile  function  F(x,y/h)  is  same  as  before.  The  value  of  Tj(x), 
however,  is  left  to  satisfy  the  other  boundary  condition  as  specified  by 
Equation  (2-42).  Completion  of  the  integration  process  from  y=h  to  y=0  would 
yield  the  improved  temperature  field. 

(5)  Revise  the  profile  function  according  to  the  improved  temperature  field.  F 
now  depends  on  both  x  and  y/h.  Repeat  the  iteration  process  from  step  (2). 
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Successful  convergence  of  Che  iteration  scheme  wouLd  be  indicated  by  a  stationary 
trend  of  Tj(x)  and  qj(x). 

Iterative  improvement  on  the  film  pressure  calculation  and  on  the  convective 
heat  transport,  presently  approximated  by  Equation  (2-38),  can  be  implemented  if 
desired.  The  temperature  field  in  the  fluid  film  would  directly  affect  the  fluid 
viscosity  in  the  hydrodynamic  analysis.  At  the  same  time  the  film  geometry  should 
also  reflect  thermal  distortions  of  the  surrounding  bodies.  With  tilting  pad 
bearings,  the  tilt  position  and  motion  should  also  be  adjusted  and  be  included  in 
the  hydrodynamic  analysis  in  a  manner  consistent  with  the  static  and  dynamic 
equilibrium  conditions  of  the  pad. 

The  above  proposed  iteration  scheme  can  be  extended  to  the  three-dimensional 
case  as  is  required  to  analyze  a  realistic  problem.  For  the  three-dimensional 
problem  Tj  and  qj  depend  upon  two  spatial  variables,  x  and  z.  Step  (2)  would  be  a 
marching  type  two-dimensional  integration.  Step  (4)  remains  a  one-dimensional 
computation. 


3.  COMPUTATION  TECHNIQUES 

For  approximately  25  years,  researchers  on  fluid  film  lubrication  have  made 
use  of  electronic  digital  computers  to  numerically  analyze  various  bearing  types. 
Theoretical  analysis  of  lubrication  can  be  formulated  as  a  two-dimensional  partial 
differential  equation  of  the  elliptic  type.  Accordingly,  the  applicable 
computation  methods  are  always  based  on  a  discretized  approximation,  which  may  be 
induced  by  either  a  finite  difference  formulation  or  a  variational  argument.  More 
recently,  the  finite-element  method  has  become  very  popular  because  it  adapts  to 
the  description  of  complicated  geometry  quite  conveniently.  As  a  computation 
algorithm,  the  finite-element  method  achieves  discretization  through  a  variational 
formalism  within  the  domain  of  a  finite,  although  small,  element. 

Interest  in  the  study  of  sliding  surface  bearings  with  a  very  small  local 
film  thickness,  and  with  the  inclusion  of  surface  roughness  and  wear  effects, 
makes  it  particularly  challenging  to  have  a  universally  accurate  computation 
algorithm  for  the  lubrication  theory.  Furthermore,  inherent  physical  features 
make  it  necessary  to  consider  elastic  compliance  of  the  bearing  parts,  thermal 
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distortion  of  the  bearing  surfaces,  and  the  temperature  dependence  of  lubricant 
viscosity.  These  issues  were  discussed  in  the  preceding  sections.  In  the  mixed 
lubrication  regime  the  film  thicknesses  are  on  the  order  of  the  surface  roughness 
of  the  mating  parts.  Emphasis  must  be  placed  upon  the  accuracy  of  the  film 
pressure  and  thickness  for  precise  considerations  of  surface  interactions. 

In  the  following  sections,  the  basis  for  implementation  of  a  highly  accurate 
algorithm  will  be  outlined.  First,  the  problem  considered  in  Appendix  A,  which 
neglects  end  leakage  and  thus  possesses  closed  form  solution  for  the  film 
pressure,  will  be  recast  in  terms  of  a  numerical  procedure  so  that  numerical 
accuracy  may  be  directly  verified  in  terms  of  this  particular  model  problem. 

Next,  a  procedure  is  provided  whereby  the  actual  metrology  record  can  be  used  to 
describe  the  surface  of  a  sector  pad  in  the  as-manufactured  state.  And  finally,  a 
high  accuracy  algorithm,  which  accounts  for  end-leakage  effects,  which  requires  a 
minimum  of  fine  mesh  treatment  to  describe  film  thickness  variation,  and  which 
permits  intermesh  interpolation  for  extrema,  is  documented. 


Finite  Difference  Computation  of  the  Infinitely  Long  Slider 

Appendix  A  treats  the  infinitely  long  slider  problem  bv  formally  integrating 
the  governing  differential  equation 


d_ 

dx 


E(H ®)  —  -  6pUE(H) 
dx 


=  0 


(3-1) 


yielding 


p(x)  =  6  pU 


f  E(H)dx  f  dx 

'  o  E(H^)  E  '  o  E(H*) 


(3-2) 


Cg  is  an  integration  constant  that  is  fixed  by  the  exit  condition 

p(b)  =  0 


(3-3) 
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Or,  in  the  case  of  a  converging-diverging  film  geometry,  if  the  rupture  condition 


E(H)  dx 
E(H 


-  EOT) 


b  Jo  EiH3) 


is  satisfied,  then 


CE  =  -E(H) 


xc  E(H)dx 
)  E(H3) 


[xc  dx  _ 
'  0  E(H3) 


together  satisfy  the  Reynolds-Swift-Stieber  condition  for  film  rupture.  To 
address  end  leakage  or  a  slider  of  finite  length,  it  is  expected  that  closed  form 
integrals  such  as  Equation  (3-2)  will  no  longer  be  available  so  that  a  numerical 
procedure  of  the  finite-difference  or  finite-element  type  will  be  needed.  As  a 
preparatory  effort,  the  infinitely  long  slider  problem  will  be  re-examined  in 
terms  of  numerical  computation  formalism,  which  subsequently  can  be  generalized  to 
deal  with  the  finite  length  slider. 

Consider  that  the  slider  width  be  divided  into  equal  intervals  such  that 

Ax  —  b/Nr 


And  the  mesh  points  are  designated  by  the  index  i,  which  represents  the  set  of 
integers  between  and  inclusive  of  (1,  N.p+1).  Note  that  x^  =  (i-lMx.  Equation 
(3-1)  will  next  be  recast  into 

<b  =  -  E(H 3)  +  6 fiUE(H) 


Consider  these  equations  between  mesh  points  (i,  i+1).  Define  the  local 
coordinate 


f  =  x-xt 

(3-10) 


Assume  that,  upon  making  Ax  suitably  small,  <t>  can  be  approximated  by  the  series 
expansion 


1 


<P  =  <PQ  +  +  ”  4  <P 2  + 


(3-11) 


To  the  same  order  of  approximation,  E(H^)  and  E(H)  can  be  similarly  expressed  as 

£(H3l  =  a0  +  4a1  +  is2o2+ ... 

(3-12) 


&lxUE(H)  =  6Q  +  +  -T62  + 


(3-13) 


Note  that  the  coefficients  of  Equation  (3-11)  are  unknown  constants  associated 
with  point  i,  whereas  those  of  Equations  (3-12)  and  (3-13)  can  be  determined  from 
their  respective  mesh  point  values.  Since  attention  is  directed  exclusively 
between  points  (i,  i+1),  the  above  series  expansions  will  be  truncated  beyond  0{f}. 


Equation  (3-8)  can  now  be  rewritten  as 


which  has  the  exact  solution 


60  ~J>0 


Setting  f  =  Ax,  one  obtains 


_(ti\ ^.(SWVlSf!) 

'  a,  /  'a,/  '  a„  / 


a„  4-  a,  /lx 


Equation  (3-9),  however,  also  requires 


=  0 


(3-15) 


(3-16) 


(3-17) 


Thus,  one  obtains 


(3-18) 

Note  that  4>Q  is  a  "universal  constant"  which  applies  to  all  i.  Equation  (3-18) 
can  thus  be  rewritten  as 


P>  +  i 


=  a  +  A  <t>n  +  B 

^1  l  0  l 


(3-19) 


(A^,  B^)  are  coefficients  which  can  be  computed  at  each  i.  At  the  leading  edge, 
i=l,  the  boundary  condition  is 


(3-20) 


Therefore, 


p2  ~  Al<p0  +  Bl 


(3-21) 


Upon  setting  i=2  in  Equation  (3-19)  and  making  use  of  Equation  (3-21),  one  obtains 


P3  —  +  Aj)^>0  +  (B2  +  Bj) 


Equations  (3-20)  through  (3-22)  can  be  given  the  generic  form 

P,=Ct<P0  +  Dt 


with  Cj=D^=0}  C2=Aj ,  1 ’  ^i+l=^i+^i ’  ^i+l— ®i+^i* 


(3-22) 


(3-23) 


(3-24a,b,c ) 


Equation  (3-24c)  can  be  implemented  repeatedly  until  i=NT  so  that  (Cnt+1,  D^t+i) 
are  regarded  as  known.  But  at  the  exit,  pNx+1  =  therefore. 


*0  =-°vr  +  1\+1 


(3-25) 


whereupon,  since  (C^,  D^)  are  now  known  at  all  i's,  Equation  (3-23)  allows  p£  to 
be  computed  also  at  all  i's. 

Now  Equation  (3-14)  allows  one  also  to  compute 


(3-26) 


noting  that  (ag,  bg)  are  known  constants  at  each  i.  Also 

(dp  \  \bo+biAx -*o\ 


(  dp  \  _  °o^  °iux 

dx  /  jv_  +  l  ao  +  ai 


(3-27) 


If  the  last  value  is  positive,  then  the  computed  set  of  p^  may  include  a  subset 
Pj,  j  being  inclusive  of  (Nc  +1,  N^.-l),  whose  values  are  negative.  One  can  then 
compute  (<}>g)j  such  that  (dp/dx)-  =  0,  i.e. 


Wj  =  <bo>j 


(3-28) 


Accordingly  (p)j  can  be  calculated  by  Equation  (3-23).  As  j  steps  up  from 
(Nc+1),  (p)j  will  be  initially  positive.  Let  j  =  Mc  +  1  be  the  first  point  where 
(p)Mc  ♦  1  is  negative,  then  the  rupture  point  had  been  established  to  be  in  the 
interval  (Mc,  Mc+^).  Whereupon  Equation  (3-14)  yields 


<P0  =(b0)Mc  +  (VmA 


(3-29) 


Substituting  into  Equation  (3-15),  with  4>^=0  and  Equation  (3-23)  used  for  Pmc»  one 
finds 


*ai\f 


<v* 


+  {bl>M 


+  D 


M 


- 


^°1  V  + 

c  c 


In 


(^) 


(3-30) 


which  is  the  characteristic  equation  for  and  can  be  solved  numerically. 
Subsequently  Equation  (3-29)  yields  d>g  and  Equation  (3-23)  can  be  used  to 
calculate  all  p^. 

To  determine  load  and  pressure  center,  one  must  seek  the  integrals 

j  pdi,  and  j  pt,di, 

Explicit  expressions  can  be  obtained  by  using  Equation  (3-15)  for  p  and 
integrating  formally.  In  a  similar  manner,  the  viscous  friction  can  be  computed. 

The  crowned  tilt-pad  bearing  analysis  of  Appendix  A  assumes  a  constant 
lubricant  viscosity  and  includes  no  coupling,  either  to  pad  distortions  caused  by 
elastic  deformations  from  load  or  to  thermally  induced  deformations. 

Consideration  of  these  topics  adds  considerable  complexity  to  the  required 
analyses . 
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Description  of  Thrust  Bearing  Film  Profile 

The  film  profile  of  a  sector  thrust  bearing  pad,  as  measured  from  a  nominal 
runner  surface,  can  be  described  as  the  sum  of  the  following  components: 

( 1 )  warpage ; 

(2)  pitch  and  roll  of  the  sector  pad; 

(3)  displacement  (axial)  and  misalignment  of  the  runner;  and 

(4)  waviness  of  the  runner 

Consider  the  sector  pad  shown  in  Figure  7.  Point  "P"  is  the  intersection  of  the 
bisecting  radius  and  the  mean  circumference.  The  four  vertices  are  marked 
(A,B,C,D).  M  is  the  midpoint  of  the  inner  arc  and  N  is  the  midpoint  of  the  outer 
arc.  E  and  F  are  midpoints  of  the  sides.  For  the  nominal  flat  surface,  all  nine 
points  are  on  the  same  plane.  Note  that  each  of  these  nine  points  is  designated 


Figure  7  -  Sector  Pad  of  a  Thrust  Bearing 


If  the  bearing  surface  is  warped,  the  warpage  can  further  be  defined  in  terms  of 
the  following  components: 


Radial  Crown  - 


Circumferential  Crown  - 


/  r  -  r  ■, 


These  can  be  respectively  measured  from  P  at  (M,N)  and  at  (E,F);  i.e. 

^  Lf  ‘  =  ^  V  ‘  —  hh 


nM  "  nP  ~~  nN  ■  nP  =  *6 
-  hp  =  hF  -  kp  —  ha 

Note  that  the  definitions  of  (hfl,  h^)  in  essence  fixes  the  metrological  setup. 

The  heights  of  the  corner  points  remain  to  be  defined.  The  following  additional 
relations  can  be  specified: 

hA+  hB  +  hc  +  hD  =  4  (ha  +  hb  +  hj 

hA'hB  +  hC  '  =  4  hd 

hA  -  hB  +  hc  +  hD  =  4  he 

hA'  hB'  hC  +  hD  ~  4  hf 

hc  is  crown  increment;  hj  is  edge  pitching;  hg  is  edge  rolling;  and  hj  is  radial 
helix.  A  lower  case  letter  is  used  as  a  subscript  to  indicate  a  type  of  surface 
deviation  from  the  nominal  plane.  Combining  these  effects. 

,,  =4(^)(-)k-*(— )*,  +  *(-)* 

warpage  \r  _r./\«/|  f  \  r  -  r  /  d  \  oc  /  « 


/  r  —  r  \ 

(rzf) 


V4  « 


2  2 


/  r  ~  rm  \  f6\ 
^  r  —  r  )  '  «  / 
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High  Accuracy  Numerical  Algorithm  for  Computing  Fluid  Film  Pressure 

The  load  capacity  of  a  fluid-film  bearing  (sliding  surface  type)  depends  on 
the  development  of  the  film  pressure.  The  film  pressure  results  from  the 
hydrodynamic  wedge  effect,  which  is  caused  by  the  spatial  variation  in  the  film 
thickness.  Operating  the  bearing  near  the  mixed  lubrication  condition  is  usually 
desirable.  The  film  pressure  is  then  concentrated  upstream  of  and  near  the 
minimum  film  thickness  location.  The  minimum  film  thickness  at  the  initiation  of 
mixed  lubrication  would  be  typically  one  order  of  magnitude  smaller  than  the 
prevailing  maximum  film  thickness.  Since  the  wedge  effect  is  inversely 
proportional  to  the  third  power  of  the  local  film  thickness,  numerical  computation 
of  the  film  pressure  requires  special  care  to  achieve  the  desired  accuracy.  The 
conventional  algorithms  of  either  the  finite  element  or  the  finite  difference 
type,  which  attempt  to  describe  the  film  pressure  by  simple  polynomial  functions, 
may  not  be  accurate  enough  unless  a  very  fine  computation  mesh  network  is  used. 
Since  a  fine  computation  network  usually  means  a  high  computation  cost,  one  often 
has  to  resort  to  use  of  a  nonuniform  mesh.  As  a  result,  the  overall  accuracy 
realizable  by  such  computation  schemes  is  unreliable  at  best.  Therefore,  a 
consistently  accurate  computation  method  for  the  fluid  film  pressure  remains  to  be 
developed.  Here,  a  new  algorithm  to  allow  high  precision  in  pressure  field 
computation  will  be  outlined. 

Strategy.  The  proposed  new  algorithm  addresses  two  somewhat  independent  issues  in 
numerical  field  computations:  the  sizing  of  the  computation  meshes  and  the 
derivation  of  the  high  accuracy  algorithm.  The  two  are  covered  in  separate 
sections . 

Background.  Christensen's^  mixed  lubrication  theory  for  longitudinal  roughness 
was  followed  for  the  derivation  of  governing  equations  in  a  form  that  is  valid  for 
both  journal  bearings  (of  either  the  full  cylindrical  or  the  partial  arc  variety) 
and  the  sector-type  thrust  bearing  pads.  Anticipating  eventual  consideration  of 


thermo-elastohydrodynamic  coupling,  the  Couette  and  the  Poiseuille  flux  components 
are  presented  in  a  form  that  allows  for  viscosity  variation  across  the  fluid  film. 
However,  since  the  mixed-lubrication  theory  has  not  yet  been  advanced  enough  to 
treat  variable  viscosity  conjunctively  with  surface  roughness,  any  symbolic 
reference  to  the  stochastic  evaluation  of  a  variable  viscosity  cannot  be  executed 
in  reality  except  for  the  limiting  case  of  smooth  surfaces. 

Following  Christensen's  lead,  variables  are  expressed  according  to  the 
stochastic  point  of  view  in  terms  of  values  of  expectation;  e.g.,  E  {p}  is  the 
expected  value  of  film  pressure  and  is  its  ensemble  average.  The  ensemble 
averaging  process,  however,  is  to  be  executed  in  the  two-dimensional  probability 
space,  with  the  roughness  distributions  of  the  two  surfaces  furnishing  independent 
stochastic  states  except  when  asperity  interference  occurs.  The  applicable 
governing  equations  are 


E{VJ  =  y2n  <»!  E{h)  +  («2  -  ~x)  E{JC) ^  -  E{Jp}  p£<p} 


—  d 


dX 


(3-31 


-  ri  ^ 

E{V\ =  -  -^77 n~E(p\ 
y  E{l/Jp}dy 


(3-32 


a  _  a  _2n  a 


- E{W)  +T  -z-E{W\  +  y  j-E{h}  =  o 
dx  *  dy  y  at 


(3-33 


Where  index  n  is  (0  or  1),  respectively,  for  the  cylindrical  journal  and  thrust 
bearings . 


(VFV,  V  )  are  dimensionless,  radius-weighted  flux  components  in  the  sliding 

x  y 

and  in  the  transverse  directions,  respectively.  (”^,  ”2^  are  rotational  speeds  of 
the  lower  and  upper  surfaces  normalized  with  a  suitable  reference  rate  w. 

(T,  7)  are  dimensionless  spatial  coordinates;  ”  is  normalized  by  the  local 
radius  r  and  is  actually  the  angular  coordinate,  and  the  transverse  coordinate  is 
normalized  by  a  nominal  radius  R.  t  is  time  multiplied  by  u>. 


"p"  is  the  film  pressure  normalized  by  Pq  u>  (R/C)2.  pQ  is  a  reference 
viscosity  coefficient.  (Jc»  Jp)  are  normalized,  viscosity-weighted  Couette  and 
Poiseuille  flow  coefficients;  they  are  defined  as  follows: 


(3-35 

p  is  the  local  viscosity  coefficient  and  is  dependent  on  the  temperature  profile. 

C  is  the  nominal  gap. 

hj  and  h2  are  heights  of  the  lower  and  upper  surfaces,  including  roughness 
profiles  on  each  surface,  as  measured  from  a  fixed  reference  surface,  h  =  h2  -  hj^ 
is  the  local  gap.  For  the  isoviscous  approximation, 

(3-36) 

For  the  present  discussion,  the  ambient  condition,  ‘p’  =  0  is  assumed  at  all  edges. 

Mesh  Sizing.  Selecting  the  proper  mesh  size  to  achieve  the  desired  computation 
accuracy  involves  a  number  of  considerations. 

1.  Structure  of  the  equation  -  The  governing  equation  is  second  order.  The 
simplest  nontrivial  field  distribution  would  be  a  parabolic  profile  that 
can  be  defined  numerically  by  no  less  than  three  points  in  each 
direction. 

2.  Field  coefficients  -  The  field  coefficients  are  mainly  associated  with 
the  film  thickness  profile  and  may  also  be  temperature  dependent  through 
Equations  (3-34)  and  (3-35).  If  crowning  of  the  surface  is  present, 
again  a  minimum  of  three  points  are  needed  to  permit  unambiguous 
numerical  definition  of  the  field  coefficients  in  the  corresponding 
direction. 


3. 


Slenderness  ratio  -  The  slenderness  ratio  for  the  cylindrical  journal 
bearing  is  L/(“R);  that  of  the  sector  thrust  pad  is  [(In  (R()/Ri))/a*.  If 
the  slenderness  ratio  is  very  large,  variation  along  the  transverse 
direction  is  minimal  except  near  either  end.  If  the  slenderness  ratio 
is  very  small,  the  short  bearing  solution^  applies,  except  near  the 
leading  and  trailing  edges.  For  a  finite  but  large  slenderness  ratio,  a 
reasonable  rule  to  follow  is  to  maintain  equidistant  mesh  spacings  in 
the  end-edge  regions  while  requiring  the  mesh  aspect  ratio  to  stay 
within  the  range  of  0.75-2.0.  For  a  small  slenderness  ratio,  mesh 
points  are  to  be  spaced  evenly  near  the  inlet  and  exit  edges.  Table  2 
illustrates  the  direct  application  of  this  rule,  assuming  three  mesh 
points  are  equally  spaced  along  the  smaller  dimension. 


TABLE  2  -  MINIMUM  RECOMMENDED  MESH  POINTS 


Slenderness 

Ratio 

Mesh 

X 

Points 

y 

Mesh  Aspect 
Central 

Ratio 

End-Edge 

1/4 

7 

3 

2.0 

1.0 

1/3 

6 

3 

2.0 

1.0 

7 

3 

1.0 

1.0 

1/2 

4 

3 

2.0 

1.0 

5 

3 

1.0 

1.0 

1 

3 

3 

1.0 

1.0 

4/3 

4 

5 

1.33 

1.0 

3/2 

4 

3 

1.0 

1.0 

4 

3 

1.5 

0.75 

2/1 

3 

4 

2.0 

1.0 

3 

5 

1.0 

1.0 

3/1 

3 

6 

O 

• 

<ni 

1.0 

3 

7 

1.0 

1.0 

*  «  is  the  pad  arc;  for  a  full  cylindrical  journal  «  =  2n 
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In  order  to  realize  good  accuracy,  the  number  of  mesh  points  in  either  direction 
probably  should  be  no  less  than  four. 

Overall  Accuracy.  The  purpose  of  an  accurate  algorithm  is  to  seek  a  numerical 
approximation  to  Equations  (3-1)  and  (3-3)  so  that  the  residue  of  the  approx¬ 
imation  vanishes  uniformly.  If  one  presumes  that  the  convergence  of  the  residue 
to  zero  can  be  achieved  by  the  limiting  process  of  an  infinitesimal  mesh  set  (Ax, 
Ay),  then  such  an  approximation  can  be  constructed  by  the  following  scheme. 


Let 


=  F(x,y)  4-  0  { residue } 
£{{£y  =  G(7,J)  4-  0  {residue} 


(3-3] 

(3-3f 


Then  F  and  G  may  be  replaced  by  their  respective  truncated  Taylor  expansions  from 
a  suitable  mesh  point  (x£,  yj);  i.e., 


dF  dF 

F(x,y)  =  Fix  y)  +  (x-x  )  —  +  iy  -  y  )  —  +  . . . 
1  J  *  dot.  J  dy 


dCr  dG 

G(x,y )  =  Gix  ,y  )  +  (x  -  x  )  —  +  (y  —  y  )  —  + 
‘  J  «  doc  J  dy. 

JJ 


Similarly,  one  can  write 


d  -  „  dH  3H 

-  E{h)  -  H(x  y  ;  +  &  ~  x  )  —  +  iy  -  y  )  —  + 

dt  1  J  1  dx  J  dy 


(3-3S 


(  3— AC 


(3-41 


Substituting  the  above  expansions  into  Equation  (3-33)  and  collecting  terms  of 
like  order,  one  obtains 


2nfjn~1Hj 


+  ...  =  0 


(3-42) 


Listing  separately  in  progressively  higher  orders,  the  approximation  requires 


dF*  dG  On 
0{1}:  — +  y"  —  +  y~nH=0 
dx  J  d y  J 


(3-43a) 


W  f  „  V  W  ()n  Oil 

0{(x  — x  )}:  —  +  ?——  +  f"^T  =0 
1  J  <9x  dv  ^  dx 


£F  nd2G  _ ndH  __,dG  9n  , 
0{(y-y)}:  ——  +  y  —  +  f  —  +  ny  —  +  2n/  // 

J  3y2  1  fyj  1  & i  J 


(3-43b) 


(3-43c ) 


Since  (dF/dxi,  dG/dyj,  H;  etc.)  are  fixed  at  the  mesh  point  (x{,  y  j ) ,  these 
expresssions  are  algebraic  constraints;  they  are  not  differential  equations. 

In  a  conventional  finite  difference  algorithm,  the  truncated  series 
expansions  are  extended  to  E{p}  such  that  Equations  (3-31)  and  (3-32)  are  directly 
reduced  to  algebraic  constraints.  In  the  proposed  approach,  Equations  (3-31)  and 
(3-32)  are  retained  as  ordinary  differential  equations  so  that  E{p}  is  not  required 
to  be  equally  smooth  as  E{TjTx}  and  E{Tpy}  Thus,  the  effects  of  variable  field 
coefficients  can  be  naturally  accommodated  in  the  locally  exact  field  of  E{p}.  In 
this  sense,  mesh  sizing  is  used  to  ensure  the  validity  of  the  approximation 
imposed  by  Equation  (3-43);  meanwhile,  within  the  domain  of  convergence  permitted 
by  the  mesh  sizes,  E{p}  is  determined  as  accurately  as  possible  through  analytical 
means . 

Polynomial  Approximations  of  Field  Coefficients.  In  the  isoviscous  case  of  smooth 
surface  problems, 

Polynomial  Approximations  of  Field  Coefficients.  In  the  isoviscous  case  of  smooth 
surface  problems, 

E{h]  =  (hi O;  E{J}  =  -  (h/O 
^  2 


vivivy1 
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E{J  }  =  —  (h/C)3\  — =  —  (WO3 
r  12  E{l/JJ  12 

(h/C)  can  usually  be  adequately  approximated  by  a  second  order  polynomial  within 
the  domain  of  two  consecutive  mesh  intervals,  even  when  the  surfaces  are  nearly  in 
contact.  Note  that,  consistent  with  a  second  order  polynomial  for  E{h},  E{Jp}.  should 
be  regarded  as  a  sixth  order  polynomial.  Such  functional  behaviors  of  the  field 
coefficients  will  be  retained.  Specifically,  one  shall  accept 

a)  E  {h}  +  <c 5  )E{J  }  =  a  +  a  (x-x  )  +  —  (x-5  )2 
i  i  i  c  o  l  t  2  1 


(3-44) 


s  =-32 


(3-45) 


Eil/JJ 


-  -  2  _  _  9 

=  C0  +  Cl^~^)+ 


(3-46) 


(aQ,  aj,  a2)  is  a  set  of  coefficients  applicable  to  the  known  data  of  E{h}  and  E  {Jc} 
f°r  *i_i  s  5  *i  +  l  and  y  =  y j  *  <b0»  bl»  b2)  is  similarly  related  to  lE{Jp}|* 
along  the  same  line  element,  (cq,  cj»  c2)  represents  the  second  order  curve  fit 
of  IE{1/Jp})'*  for  x  =  xi,  y j_ i  ^  9  s  y j+i  - 


Locally  Exact  Ordinary  Differential  Equation  Solutions.  Consider  the  2x2  mesh 
cluster  around  point  (x^,  yj)  shown  for  the  general  case  in  Figure '8.  Uniform 
mesh  intervals  are  not  assumed.  Based  on  the  known  numerical  data  of  the  field 
coefficients,  (ag,  a^,  a2;  bg,  b^,  b2;  eg,  c^,  c2)  can  be  obtained  by  second  order 
curve  fitting. 


Substituting  Equation  (3-39)  for  E{q7x}  in  Equation  (3-31)  and  making  use  of 
Equations  (3-44)  and  (3-45),  one  obtains 


dE$} 

dx 


+  . . . 


IVV*-*V  +  aV*-”,*2!3 
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(3-47) 


Similarly,  substituting  Equation  (3-40)  for  E{qJy}  in  Equation  (3-32)  and  making 
of  Equation  (3-46),  one  obtains 


dE{p} 

$ 


dG 

-G  ~  — —  (y-y)  - 
ij  37  ‘ 

J  j 


1 yJ  +  (y-yJ)\nk0+cl{y-yj)+$c2(y-y)2]:i 


(3 


Ay(0 


use 
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Figure  8  -  Schematic  of  a  Nonuniform  Mesh  Cluster 


Closed  form  integration  of  Equation  (3-47)  along  y  =  y-  from  x  =  x^  yields 

/  dF.  .  ^ 

mr  IA0(»-A0(ri  Kjf.0-F,  j(  +  1A,(B  -  Aia,ll(jfa1-  j  +  .  . . 


(3-49) 


where , 


i0(x)  = 


Li(5)=  f 


2  o 

6.4-  b(x-x  )  +  —  (x-x  ) 
01  1  2  1 


(x  —  x  )dx 

L 


2  2 

6.4-  6,(1— x  )  4-  —  (x—x  ) 

0  t  i  2  ‘ 


(3-50a) 


(3-50b) 


The  detailed  expressions  for  (Ag,  )  are  given  in  Appendix  B.  Similarly,  closed 
forms  integration  of  Equation  (3-48)  along  x  =  x£  from  y  =  yj  results  in 


E{p}~  =  -lB0(y)  -  B0(Jj)\Gi  ;  -  (B,(y)  -  Bft) i  - 


(3-51) 


where , 


BQ(y)  = 


(y)n{c0+  ctf-y)  4-  c2(J-y;)2/2}3 


r  (y-y,)dy 

1  *  (y)"{c0+  ct(y  -  y)  4-  c2(y  -  „7  )2/2}3 


( 3-52a ) 


(3-52b) 


Details  of  (Bo,  B],)  are  given  in  Appendix  C. 

The  Central  Algorithm.  Truncating  the  expansion  of  E{q»x}  at  the  first  order  term, 
(3-49)  can  be  regarded  as  containing  two  undetermined  coefficients  (Fi,j  and 
3Fifj/8xi).  The  data  of  E{p}  at  (i-1,  j)  and  (i  +  1,  j)  can  be  utilized  so  that 


1  -1  O'!  \E\pi.i  Jf  I  r Aofxt-i)  •  A0(xi)  Ai(xi.i)  -  Aj(xi)  "I  \yfaa0-Fitj 


0  -1  1  J  I  E{pi  +  l  ,j}  I  lAofxt  +  i )  -  A0(xi)  Atixi  +  i)  -  Ai(xi) J  ly?"  at  -  dFiJdxi 


The  matrix  equation  can  be  inverted  to  yield 

\A\  l yfnao-FiJ  =  {Ai(*t+i )  'At(*i)JE{Pi.iJ-fAi(sn.i>-AI(ii.i)E{piJ 

-  [Ai  (*i_t)  -  Aj  xj}  E  {pi+i  j} 


|a|  {yfnat -dFJdx}  =  -[A0(xi+l)  -AgfxJ] Efp^j}  +  [A0(xi+I)  -Ajx^)] E{pJ 

+  [A0(xi  t)  -A0(x)]E{pi+IJ 


where 


A  - 


Aj x.  , )  -  AAx)  AAx  ,)-AAx) 
0  c  —  1  Oi  1  *  —  1  li 

AAx..)-AAx)  AAx)~AAx) 
Oi+I  0  t  lj+1  1  i 


Similarly,  making  use  of  the  data  of  E{p}  at  (i,  j-1)  and  (i,  j+1)  would  lead 
to 


Wj-i'-W  'V 

BJy.  .  ,)  -  BAy  A  BAy.  ,  ,)  -  BAy  A 


To  simplify  symbolics  further,  define 


(3-60a,b) 

Combining  (dF^j/djcj,  dGijj/dyj)  respectively  from  Equations  (3-55)  and  (3-58)  into 
Equation  (3-43a),  using  the  new  notation,  it  is  found  that 

[E{p  ,  .} 

+  -Ao<  +  >+A<x->-Bo<  +>  +  B<x  -t  -AoJ  w.  }  +  = 

k.) 

(3-61) 

This  is  the  high-accuracy  central  difference  algorithm.  It  is  to  the  variable 
viscosity  lubrication  theory  for  surfaces  with  longitudinal  roughness  as  the 
relaxation  formula  is  to  the  classical  potential  theory.  Its  inherent  accuracy  is 
not  affected  by  the  presence  of  localized  small  film  thickness  provided  the 
surface  profile  is  accurately  described  by  a  second  order  polynomial. 

Equation  (3-61)  applies  to  all  (i,  j).  If  one  chooses  a  fixed  i,  the  total 
system  encompassing  all  j  can  be  written  as 


[(/«]{p,./}  +  [V.Kp,}  +  =  {<?.} 


(3-62) 


{pi}is  a  column  vector  with  elements  of  E{p^  j}  (fixed  i,  ail  j).  (U£)  is  a 

diagonal  matrix  with  the  element  Aq(+).-  (V^)  is  a  tri-diagonal  matrix  with  the 
elements  ( Bo ( ♦ ) »  "  A0(+)  +  A0(-)  -B0(+)  +  B0(-)»  "  B0(-)*  ^Wi^  another  diagonal 
matrix  with  the  elements  -Aq(_j.  {Q^}  is  a  column  vector  containing  the  elements 

The  rank  of  the  system  represented  by  Equation  (3-62)  is  the  number  of  mesh 
points  in  the  y-direction.  Inversion  of  this  system  is  the  bulk  of  the 
computation  effort.  For  a  bearing  with  a  large  slenderness  ratio,  there  are  fewer 
mesh  points  in  the  x-direction.  Therefore  it  would  be  more  advantageous  to  treat 


the  alternative  representation  (fixed  j,  all  i),  which  can  be  written  in  the  form 

+  W/Pj  +  i*  =  IQ} 


(3-6 


Determination  of  the  Rupture  Boundary.  The  presence  of  film  rupture  is  indicated 
by 

£{p}  -  pc 

( 3-6' 


whereupon  the  rupture  boundary  is  located  by  the  condition 

E{p\  =  P„ 


(3-6: 


From  the  standpoint  of  numerical  computation,  where  E{p}  is  represented  by  mesh 
arrays,  a  sufficient  condition  for  the  inequality  (3-64)  is 


(3-6i 


Equation  (3-65)  can  be  solved  iteratively,  in  accordance  with  Equation  (3-49)  or 
(3-51),  to  establish  the  intermesh  location  of  the  rupture  boundary.  Once  the 
rupture  boundary  is  located,  then  one  may  seek  to  impose  the  Reynolds-Swif t- 
Stieber  condition  as  the  criterion  for  subsequent  global  iterations. 

If  the  inequality  (3-64)  occurs  in  the  interior  of  four  surrounding  mesh 
points  then  a  more  elaborate  scheme  is  needed  to  determine  the  rupture  boundary. 
Since  this  situation  is  not  likely  to  be  associated  with  the  thin-film  operating 
condition  of  sliding  surface  bearings,  it  will  not  be  given  further  attention. 


4.  SURFACE  ROUGHNESS  EFFECTS 

When  a  sliding  surface  bearing  operates  with  a  minimum  film  thickness,  which 
is  of  a  similar  order  as  the  surface  roughness,  the  classical  lubrication  theory 
is  no  longer  adequate.  Allowance  must  now  be  made  for  asperity  interactions  — 
load  is  shared  between  hydrodynamic  pressure  and  asperity  contacts.  At  asperity 
contacts,  friction  behaves  according  to  the  principles  of  boundary  lubrication, 
and  attendant  heating  and  wear  inevitably  lead  to  alterations  in  the  surface 
texture.  The  role  of  surface  topography  in  tribological  functions  cannot  be 


expected  to  be  uniquely  associated  with  a  single  decriptive  parameter,  e.g., 

C.L.A.  (center  line  average)  or  R.M.S.  (root-mean  square)  roughness.  Researchers 
are  beginning  to  bring  together  topics  that  had  been  regarded  as  isolated  issues 
to  examine  the  detailed  mechanisms  of  these  complicated  interactions. 

Christensen's  theory  of  mixed  lubrication  for  longitudinal  roughness  has  been 
selected  as  the  analytical  point  of  view,  since  the  longitudinal  orientation  of 
topographical  features  seems  to  be  a  reasonable  model  for  machined  and  worn-in 
surfaces. 

In  Appendix  A,  as  in  Christensen's  "A  Theory  of  Mixed  Lubrication, the 
model  concentrates  on  the  behavior  of  the  bearing  with  longitudinal  roughness  on 
one  surface  only.  Christensen  points  out  that  if  both  the  surfaces  in  a  bearing 
are  rough,  then  the  probability  density  function  used  will  be  considered  as  a 
description  of  the  combined  roughness.  The  roughness  distribution  functions  used 
by  Christensen  were  symmetrical  and  approximated  the  Normal  or  Gaussian 
distribution.^’^  There  remain  areas  for  theoretical  improvement  such  that  more 
realistic  characterization  of  the  interacting  surfaces  may  be  considered.  The 
subsequent  discussions  will  deal  with  two  specific  aspects  of  conceptual 
generalizations  that  are  related  to  the  stochastic  description  of  the  topographs. 
First,  the  roughness  deviations  on  each  surface  are  recognized  as  an  independent 
set  of  statistical  data;  consequently  a  two-dimensional  representation  of  the 
stochastic  model  would  have  to  be  used.  Second,  a  machined  surface  or  one  that 
had  experienced  prior  wear  would  have  fewer  peaks  than  valleys;  thus  the  roughness 
deviations  should  have  skewed  distribution  functions.  If  these  improvements  were 
incorporated  into  mixed  lubrication  theory,  all  stochastic  operations  used  in  the 
theory  would  have  to  be  consistently  revised. 


Mixed  Lubrication  of  Sliding  Surface  Bearings  with  Long’tudinal  Roughness 
•  3 

Background .  Christensen  has  proposed  stochastic  models  for  hydrodynamic 
lubrication  of  rough  surfaces  and  has  applied  this  technique  to  study  the  mixed 
lubrication  problem  using  the  one-dimensional  flat  slider  for  illustration.^ 
Christensen's  theoretical  models  include  striated  roughness  of  both  transverse  and 
longitudinal  orientations  as  well  as  a  uniform  isotropic  roughness.^  He 
considered  the  hydrodynamic  aspect  using  the  thin  film  model,  which  is  strictly 
applicable  only  with  roughness  on  the  stationary  surface. 

Q 

Elrod  considered  striated  roughness,  again  using  the  thin  film  model,  from 
the  deterministic  point  of  view.  The  latter  approach  is  a  direct  extension  of  the 
lubrication  theory  for  bearings  with  shallow,  narrow  grooves  on  one  surface 
only.^’^®  In  the  most  general  form,  both  surfaces  may  be  rough;  however,  their 
respective  striations  must  be  parallel.^  Using  a  multiple-scale  expansion 
method, 12  it  was  shown  that  the  flow  across  the  striations  would  not  follow  the 
rapid  gap  fluctuations  associated  with  surface  roughness.  If  the  surface  motion 
has  a  component  that  is  transverse  to  the  striations,  then  the  moving  roughness 
would  contribute  to  an  additional  squeeze-film  effect,  which  has  a  nonvanishing 
temporal  average.  Patir  and  Cheng,  concerned  mainly  with  transverse  striations, 
arrived  at  similar  conclusions.  More  recent  efforts  on  lubrication  analysis  for 
rough  surfaces  are  reported  in  References  14,15  and  16.  In  Reference  17,  Elrod 
gave  a  comprehensive  review  of  theoretical  ideas  on  roughness  effects,  including 
deterministic  considerations  of  two-dimensional  Reynolds  roughness  (thin  film), 
Stokes  roughness  (short  wavelength),  as  well  as  molecular  slip  (very  thin  gas 
films) . 

The  purpose  of  this  section  is  to  consolidate  relevant  approaches  to  form  the 
basis  for  studying  mixed  lubrication  on  both  stationary  and  moving  surfaces  to 
allow  for  longitudinal  roughness.  Specifically,  equations  suitable  for  the  study 
of  cylindrical  journal  and  sector  type  thrust  bearings  will  be  developed. 


Description  of  Surfaces.  Roughness  Functions,  and  Film  Thickness.  Subscripts  1 
and  2  will  be  used  to  designate  the  lower  and  upper  bearing  surfaces, 
respectively.  Each  point  on  the  surface  is  measured  as  the  height  from  an 
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appropriate  nominal  surface;  and  each  surface  is  described  as  the  sum  of  a 
geometrical  part  (which  is  "smooth")  and  a  topographical  part  (which  is  "rough"). 
Thus 

hj  =Ht  (x/b,y/L,t)  +8J(y/\I) 
h2  =  H2  (x/b,  y/L,t)  ■  82  (yl \2) 

(4-la,  b) 

(x,y)  are  local  Cartesian  coordinates,  x  being  aimed  in  the  sliding  direction;  t 
is  the  time.  8^  and  §2  are  topographical  descriptions  of  the  roughness  features. 
Their  sign  conventions  depict  elevation  above  the  respective  mean  surfaces  as 
positive  values.  (Ai,  are  roughness  scales  for  the  respective  surfaces.  It 
will  be  assumed  that 


(tvt2)  =  (\t,\2)/L  <  <  1 


(a,,  a2) 


are  rms  values  of 


a 


2 

l 


I 


O  1  A  * 

y0-hAy  Ay 


J  y  -  $  A  y  Ay 
o 


It  is  assumed  that 


(\,,\2)  <  <  A y  <  <  L 


(4-2) 


(4-3a,  b) 


(4-4) 


(oj,o2)  can  be  dependent  upon  (x,  y),  if  desired,  to  allow  for  spatial  variation  of 
roughness  features. 

Separation  of  the  surfaces  is 


h  -  h2-ht 
=  H-  S.-S, 

(4-5 ) 
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where 


h  =  h2-hi 


(4-6  ) 


is  the  usual  film  thickness  upon  neglecting  roughness.  5j  and  82  may  be  regarded 
as  definable  through  their  respective  profilometry  data.  However,  the  coordinates 
(y/Aj,  y/A2)  can  not  be  related  to  each  other.  In  this  sense  (8j,  S2)  should  be 
regarded  as  stochastically  independent  of  each  other;  that  is,  the  probability  for 
8j  to  assume  any  particular  value  is  independent  of  the  particular  value  of  82. 

This  point  of  view  is  crucially  important  to  the  development  of  a  consistent 
stochastic  analysis. 


Stochastic  Computations.  The  stochastic  descriptions  of  the  topographic  functions 
(Sj.Sg)  will  be  assumed  to  be  known  in  terms  of  their  respective  surface  height 
histograms  P^  (8j)  and  ?2  (S2).  Pjd  8j  can  be  regarded  as  the  infinitesimal  area 
fraction  on  the  lower  surface  that  has  the  (roughness)  height  of  8,.  (P^,  P2)  are 

sometimes  called  probability  density  functions.  By  the  nature  of  this  definition, 
the  following  identities  hold: 

e  ao  r  ao 

V*i  =  l;  p2dS2  =  1 

J  —00  J  —CD 

(4-7a,  b) 

r  ao  r  ao 

I  W^  =  0;  P2S2dS2=  0 

J  —CD  J  —CD 

(4-8a,  b) 

®?=  f  PAdSv  °l=\  P2S\dS2 

J  _ ao  J  _ ao 

(4-9a,  b) 


The  upper  and  lower  limits  of  all  real  roughness  profiles  are  finite;  therefore, 
the  integration  limits  of  Equations  (4-7)  to  (4-9)  should  be  replaced  accordingly. 
(8,,82)  are  independent  variables  of  the  phase  space  in  the  stochastic  analysis  of 
a  problem  involving  two  rough  surfaces.  The  expected  value  of  a  function  f  that 
is  stochastically  dependent  on  (8,,82)  is  computed  by  ensemble  averaging 
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(integrating  with  probability  weighting)  in  the  two-dimensional  stochastic  phase 
space.  Formally  this  is  written  as 


E(f)  =  /  J(P,  P2)  f(d&y  d52) 


(4-10) 


The  applicable  domain  excludes  the  condition  8,  ♦  82  >  H,  since  H  is  the  separation 
of  the  mean  surfaces.  The  area  fraction  of  asperity  contact  is 

Ac  =  ff  ( P,P2)  (dSt  dS2)  > 0 


8,  +  S2  >  H 


(4-11) 


General  Equations.  The  isoviscous,  laminar  film  theory  of  lubrication  will  be 
accepted  for  the  present  treatment.  In  vector  notation,  this  is 


diu  V  +  —  =0 
dt 


(4-12) 


(4-13) 


V  is  the  volume  flux  vector.  V  is  the  sum  of  the  sliding  vectors.  This  form  can 
be  used  for  a  counter-rotating  bearing.  For  cylindrical  journal  and  sector  type 
thrust  bearings,  these  equations  can  be  rewritten  in  the  component  form  with  a 


common  representation. 


x  =  rQ;  r  =  R  (y!R)n 


(4-14a,  b) 


h  (toj  +  r 


_1_ 

12 fir  dy' 


if/  4  dp 

y  12  fi  dy 


(4-15a,  b) 


3E 


n=(0,l)  respectively  for  the  cylindrical  journal  and  thrust  bearings. 
For  a  dimensionless  formulation,  the  spatial  coordinates  are 


x  =  9  =  x/r  ~y  =  y/R 


(4-17a,  b) 

A  nominal  clearance  C  is  used  to  scale  the  composite  film  thickness  while  the  two 
roughness  profiles  are  respectively  scaled  by  their  rms  roughness  values  (oj,o2): 

h  =  hJC  =  H  - “j  8j  -0*2 82 


(4-18) 


where , 

H^(H2Hj)IC 

(EJ,TT2)  =  (Bj/Oj.e^);  (oj.o^)  =  (oj/C.Og/C) 


The  rotational  speeds  are  scaled  by  a  nominal  rotational  rate 

Cw1(&>2)  =  (tUj/w,  u>2/«) 


The  flux  components  are  scaled  by  the  nominal  Couette  flux 

(VxfVy)  =2(WX,  'Py)/(0)fiO 


(4-19a,  b,  c) 


(4-20) 


(4-21) 


The  film  pressure,  as  measured  above  ambient  pressure,  is  scaled  by 


Equations  (4-15,  4-13)  can  thus  be  rewritten  in  the  dimensionless  form: 


V  =  Go.  +aJhy"-  —  V  =  -T3  — 
x  1  *  J  ST  y  gy- 


(4-23a,  b) 


1  3  1  3  _  .  84 


*  u  »  Ur* 

+ - (p"?  )  +  —  =  o 

3T  T"  if  y  dt 


(4-24) 


The  natural  dimensionless  time  is 


t  =  u>t/2 


(4-25) 


Multiple-Scale  Analysis  of  Striated  Reynolds  Roughness.  After  Elrod, ^  one 
recognizes  that  the  roughness  function  is  more  appropriately  described  in  terms  of 
a  local  coordinate  of  a  fine  structure;  i.e., 


«i-v(rzK-v 


1  2 ' 


(4-26) 


Therefore,  partial  spatial  derivative  with  respect  to  the  transverse  coordinate 
must  also  make  allowance  for  these  fine  scales,  which  play  the  role  of 
amplification  in  the  differentiation  process.  In  other  words 


a_ 

W 


_  £ 

<r 


+  (R\  £_  I  +(R\ 

“<^<2  It“<2  '  ^2 


^2  *•  ”<| 


(4-27) 


(R/Aj,  R/A2)  can  be  regarded  as  roughness  amplification  factors  associated  with 
transverse  (relative  to  the  direction  of  striation)  spatial  differentiation. 


One  may  also  reasonably  postulate  that  longitudinal  spatial  differentiation 
and  temporal  differentiation  do  not  introduce  roughness  amplifications.  A  simple 
idea  is  used  in  the  narrow  groove  lubrication  theory  and  is  known  to  be  valid  in 
all  areas  remote  from  the  edges.  Equation  (4-24)  can  thus  be  rewritten,  treating 
7lt  £j,  <2  as  independent  variables,  and  be  regrouped  according  to  various  degrees 
of  roughness  amplification: 


(4-28 


V  +  — 

y  -n 


X 


+  -  (yn  V  ) 
dy  y 


+ 


an 

-=.=  0 

dt 


This  expression  suggests  the  usefulness  of  constructing  the  solution  in  terms  of  a 
double  power-series  expansion  with  respect  to  the  two  smallness  parameters 
(A,/R,  \2/R). 


Neglecting  o{Aj/R,  A2/r},  since  Equation  (4-28)  is  to  be  valid  for  the 
independent  limits  (Aj-»0,  A2-»0),  one  concludes  that 


dV 


*1 


=  0; 


dV 


*2 


=  0 


(4-29 


These  expressions  are  the  essence  of  the  theory  of  striated  Reynolds'  roughness. 
The  flux  component  transverse  to  the  roughness  striations  does  not  respond  to  the 
fine-scale  fluctuations  of  the  surface  topography.  This  general  conclusion  is 
independent  of  the  orientation  of  the  striations  relative  to  the  sliding 
direction,  but  does  require  the  striations  on  one  surface  be  parallel  to  those  on 
the  other. 


Stochastic  Analysis  of  Film  Flux  and  Film  Pressure.  The  ultimate  interest  of 
lubrication  analysis  is  concerned  with  the  spatial  integrals  of  film  flux  and  film 
pressure.  In  the  statistical  studies,  these  are  really  ergodic  averages  as 
opposed  to  ensemble  averages.  Thus  for  a  specific  rough  surface,  spatial 
integration  is  the  real  objective,  and  stochastic  analysis  is  conceptually 
irrelevant.  If,  however,  the  surface  topography  can  be  defined  only  in 
statistical  terms,  then  the  results  are  valid  for  all  surfaces  that  are 
statistically  similar  to  the  sample  topography.  In  this  sense,  the  statistical 
point  of  view  begins  to  have  meaning.  However,  the  use  of  ensemble  averaging 
instead  of  spatial  integration  is  mainly  a  matter  of  convenience.  When  both 
surfaces  are  rough,  because  the  relative  location  of  the  two  roughness  profiles  is 
not  certain,  the  actual  composite  roughness  cannot  be  uniquely  constructed  from 
the  individual  roughness  profiles  in  a  deterministic  manner.  Accordingly,  the 
stochastic  point  of  view  becomes  a  matter  of  necessity.  In  the  following, 


stochastic  analysis  will  always  refer  to  ensemble  averaging  (expectation)  and  the 
deviation. 


Symbolicly  every  function  will  be  so  treated;  i.e., 

f=E(f)+f 


(4-30) 


The  expectation  operation  has  already  been  defined  by  Equation  (4-10),  in  which 
the  roughness  heights  (8j,82)  were  used  to  define  the  probability  density  functions 
(Pj,  P2).  Consistent  with  the  dimensionless  notations  adopted  above,  the 
expectations  operation  on  a  dimensionless  variable  will  be  defined  in  terms  of  the 
dimensionless  probability  density  functions. 

/^(Sj)  =  Oj  Pj  (S^OjSj); 

P2(82)  3  o2  P 2  (82  =  o282)  (4-31a,  b) 

# 

Clearly,  expectations  of  the  film  fluxes  and  the  film  pressure  are  of  particular 
interest.  Equation  (4-29)  is  equivalent  to  a  null  condition  of  the  deviation 
component  of  'Pv,  i.e., 

V  =  V  -  E(W  )  -  0 

y  y  y 

(4-32) 


Equation  (4-23b)  now  can  be  rewritten  as 


®  V 


(4-33) 


Taking  the  ensemble  average, 

(4-34) 


or, 
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f 


Equation  (4-33)  can  be  differentiated,  with  respect  to  x,  yielding 


l(&)=  _ 
ar\a7/  arvT3/ 


(4-3t 


The  order  of  differentiation  of  the  left-hand  side  can  be  written  in  reverse; 
therefore , 


±(d2)=  _  ±(^y) 

Sy'Kdx)  dx  \  T3  / 


(4-37 


The  last  relation  means,  neglecting  0{X^/R,  A2/R}’  that  dp/dx  also  has  a  null 
deviation  component;  i.e. 

^=0 

dx 

(4-3E 

Consequently,  from  Equation  (4-23a),  the  ensemble  average  of  Vx  can  be  derived  as 

_  dE(p) 

ryttrt  \  /*  \  Wl  r:»/  L  \  r/T«i\  r 


E(V  )  =  (S>  +  S7Jy”  fifli  -  —  E(P)  -±- 

*  1  2  J  —n  ,tF 


(4-39 


Finally,  taking  the  ensemble  average  of  Equation  (4-24)  term  by  term,  one 


obtains 


r 


(4-40 ) 


Computation  of  Ensemble  Averages.  To  solve  (4-40)  for  E  Cp),  ensemble  averages  of 
(F,  ,  1/h^)  are  needed.  The  procedures  of  these  computations  depend  on  whether 

or  not  asperity  contacts,  would  take  place.  Let  (-aj,  b^)  be,  respectively,  the 
(lower,  upper)  limits  of  the  roughness  profile  of  the  lower  surface,  and  {-&2*  ^2^ 
be  those  of  the  upper  surface.  Then  the  condition  for  asperity  contact  is 

"ojbj  +  <T2b2  >  H 

(4-45) 
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This  inequality  marks  the  boundary  of  a  forbidden  domain  in  the  two-dimensional 
phase  space;  this  domain  must  be  excluded  from  all  ensemble  average  calculations. 
Figure  9  illustrates  these  ideas. 


Figure  9  -  Applicable  Domain  for  Ensemble  Averaging 


It  is  convenient  to  define  the  unrestricted  expectation  as  the  ensemble 
integration  over  the  entire  rectangular  domain: 


r 6!  _  _  r  *2  _  _ 

£</)=  PtO[l  P2fdS2 

J  -a,  J  -an 


(4-46 ) 


Then,  if  the  inequality  (4-45)  is  satisfied,  also  define  the  expectation 
deficiency  as 


4£v>  =  f,JT. -Tfl  V<«2 


(4-47) 


The  net  expectation  is 


E(f)  =  Ea(p.AE(fi 


(4-48) 


For  the  various  functions  of  interest,  the  unrestricted  expectations  are 


E0(h)  =  H 


Eo(P)  =  H3  +  3  (Oj2  +  o^TT-o;3  Eo  (5^3)  -F23  Eo  (823) 


(4-49a,  b) 


E0(l/F3)  is  much  more  complicated  and  will  be  given  attention  separately. 

Expectation  deficiency  functions  of  CK’j’K'3)  require  computation  of  the  type 


T'P’dfi  n  =  0,1, 2, 3 


through  appropriate  ranges  as  defined  in  Equation  (4-47).  In  particular 


(4-50) 


dBU)=  j(7r__  l  .  r,dd,  —  t\  rndd 


2b2,  1-1  J  (ff-Ojffj.  ‘  2  2 


(4-5: 


yields  the  area  fraction  of  contact. 


Hydrodynamic  Shear.  The  hydrodynamic  shear  stress  in  the  sliding  direction  is  of 
interest  because  it  contributes  to  the  drive  torque.  Depending  on  the  particular 
circumstance,  the  shear  stress  on  either  the  upper  or  the  lower  surface  is  of 
interest;  they  are 


M(w2-o}l)r  +  hdp 
h  ~  2  dx 


(4-5: 


The  upper  sign  is  to  be  used  for  the  shear  stress  on  the  upper  surface;  and  the 
lower  sign,  for  the  shear  stress  on  the  lower  surface. 

Scaling  with  pwR/C 


T=  5X^7*  ±3 

h  Zn  c>7 


(4-5: 


Carrying  out  ensemble  averaging, 


—  _  E(h )  dp 

E( r)  -  (o)  —  co  )y  E(l/h)  ±  3  -  — 

1  1  yn  dx 


(4-5< 


EQ/TT)  is  the  additional  expression  that  did  not  appear  in  the  analysis  of  film 


pressure . 


Stochastic  Computations  for  Mixed  Lubrication  with  Skewed  Roughness  Distribution 
Function 

The  theory  of  mixed  lubrication  —  for  example,  the  original  treatment  by 
Christensen^  —  was  developed  with  a  symmetrical  roughness  distribution  function. 
The  bearing  roughnesses  in  practical  applications  are  not  likely  to  be 
symmetrically  distributed,  rather,  they  are  skewed.  Most  likely,  there  are  more 
valleys  of  a  certain  magnitude  than  peaks  of  the  same  magnitude.  The  necessary 
analysis  to  accommodate  an  unsymmetrical  roughness  distribution  function  is  given 
below. 
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Ensemble  averaging  (as  the  expectation)  of  a  stochastic  function  that  depends 
on  two  random  variables,  f  (8^,  82),  is  defined  as  the  probability-weighted 
integral  in  the  two-dimensional  phase  space;  namely. 


E(f)  =  ff[P,  (Sj)P2(S2)]f(dS,  dS2) 

With  a  finite  separation,  H,  the  domain  of  integration  should  exclude  the 
condition  of  interference,  which  is 

s,+s2>h 


(4-60 


(4-61 


Consequently,  the  ensemble  average  can  be  treated  as  the  sum  of  the  unrestricted 
expectation  and  an  expectation  deficiency; 


E(0  =  EJf)  -AE(P 


(4-62 


where , 


'.#=  f_J 


AE(P  =  //  Pj  P2fdSldS2 


S  j  +  S2  >  H 


(4-63 


(4-64 


The  probability  density  functions  (Pp  P2)  should  be  based  on  actual 
topographical  records  of  the  surfaces.  However,  to  facilitate  the  computation 
procedure,  an  empirical  function,  which  shares  the  same  overall  characteristics  of 
the  actual  function,  may  be  used.  Christensen  proposed  a  symmetrical  sixth  order 
polynomial,  which  contains  a  single  free  parameter  that  can  be  directly  associated 
with  the  rms  roughness,  o,  defined  by 

o2  =  IPS2dS 


L, 
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I'.'- 


y 


B 


as 


Christensen's  empirical  formula  is 

P  = 


32(3o)' 


[(3o)2  -  <?2]3 


P  (S2  >  (3a)2)  =  0  forS 2  <  (3a)2 


(4-66) 


Note  that  this  function  has  vanishing  derivatives  up  to  the  second  order  at  the 
extreme  points  (8=  ±3o).  Furthermore,  the  distribution  of  asperity  height  is 
symmetrical,  so  that  there  are  as  many  protrusions  as  there  are  depressions.  For 
real  surfaces,  however,  due  to  the  inherent  features  of  machining  processes  and 
possible  run-in  effects,  there  always  are  fewer  protrusions  than  depression. 
Therefore,  the  centroid  of  the  distribution  function  should  be  skewed  to  the 
negative  side  (depressions).  It  should  be  possible  to  modify  Equation  (4-66) 
slightly,  so  that  a  second  parameter  becomes  available  to  describe  skewness  in  the 
topography.  A  seventh-order  polynomial  will  be  examined  as  a  proposed  basis  for 
the  computation  of  ensemble  averages  of  interest. 


Basic  Properties  of  the  Distribution  Function.  The  distribution  function  is 
always  positive,  and  must  satisfy  two  contraints  due  to  the  very  nature  of  its 

definition.  The  first  is  the  convention  of  normalization  in  probability  theory: 

+  00 

PdS=  1 


(4-67) 


The  second  is  the  condition  of  null-biasing,  which  must  be  satisfied  because  the 
stochastic  function  itself  must  have  a  null  value  for  its  ensemble  average,  i.e., 

PSdS= 0 


(4-68) 
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The  rms  deviation,  sometimes  called  the  standard  deviation,  o,  was  previously 
given  by  Equation  (4-65).  These  are  moment  integrals  of  P(8)  of  order  (0,  1,  2) 
for  Equations  (4-67,  68,  and  65),  respectively.  One  way  to  examine  the  skewness 
of  the  distribution  is  to  consider  the  third  order  moment  integral  as  given  by 


<-r 


PS3  dS 


(4-69 


Note  that  for  a  symmetrical  probability  density  function  such  as  Equation  (4-  66), 
and  all  odd  orders,  the  moment  integrals  vanish  identically. 


Derivation  of  a  New  Probability  Distribution  Function.  Instead  of  using  Equation 
(4-66)  as  the  probability  distribution  function,  its  basic  form  is  adopted  as  a 
generating  function. 

G(8)  =  (c2  -  82)3  for  82  <  c2 


G(8)  =  0for82  >c2 


The  indefinite  moment  integral  of  G  is  defined  as 


(4-70 


c6  3  c4  S2  3  c2S*  S 6 


„  (  c°  3c  8  3 ( 

r  (S)=  GS  dS  =  (  —  -  -  +  — 

n  '  n+ 1  rc  +  3  n 


+  5  n  +  7 


(4-71 


One  may  postulate  that 


P(8)  =  [A+B  (8  +  80)|  (c2  -  (8  +  S0)2]3  for  (8  +  8Q)2  s  c2 


(4-72 


For  Equation  (4-67)  to  be  satisfied 


l  =2  A  I  (c)  = 


32  Ac' 


T 


Accordingly,  the  various  statistical  parameters,  all  normalized  with  o,  are 


Half  maximum  Height  “  =  (c/o) 

Skew-Bias  Parameter  8Q  =  (8 Q/o) 

Root  Mean  Third  Moment  9g  =  (9g/o) 

The  corresponding  probability  distribution  function  is 

7(8)  =  o 7(8) 


(4-79a,b,c) 


(4-80) 


Recasting  Equation  (4-75)  in  the  new  notation,  one  can  write 

2  3 

p(3)  =  (  ??)(=)  l+f i)T(T+7>  i-[(7+7yr  for-7-J  <7 <7 -7 

\  32  /  V  c  /  l  \“/°  »  I  °  \  °  ° 


Equations  (4-81)  and  (4-77)  can  be  rewritten  as 


7=  3  V 1  +  8* 


(4-81) 


5’=(nHK2-3) 


(4-82) 


(4-83) 


These  expressions  provide  an  easy-to-use  single  parameter  (8Q)  for  the  statistical 
characterization  of  skewed  topography. 

From  the  specific  profilometry  data,  the  value  of  9g  can  be  calculated.  Thus, 
according  to  Equation  (4-83),  the  skew-bias  parameter  can  be  determined  as  the 
root  of 


However,  because  this  solution  is  not  unique,  the  appropriate  branch  must  be 
identified.  Equation  (4-84)  lacks  a  unique  root  because  8g  has  a  stationary  point 
at  8q2  =  1/8.  Fortunately,  in  the  range 


o  «s  8 ‘s  - 
0  8 


(4-85) 


(0,  80)  nave  a  one-to-one  correspondence.  In  fact  the  upper  bound  of  802  makes 
certain  that  p  is  non-negative  in  the  full  range  of  8.  This  constraint  is 
obviously  due  to  the  relatively  simple  functional  form  of  Equation  (4-81),  and 
limits  its  applicability  to 

l-,l  V? 


(4-86) 

If  a  larger  numerical  value  of  0g should  be  indicated  by  actual  profilometry  data, 
then  the  upper-bound  condition,  (8Q2=l/8),  would  have  to  be  used  in  lieu  of  solving 
Equation  (4-84)  for  80. 

If  the  last  inequality  is  satisfied,  the  root  of  Equation  (4-83)  can  be  rapidly 
found  by  the  Newton-iteration  method.  Given  an  initial  estimate  8  o,  the  residue 
R'  =  R  (8’0)  *  0  can  be  readily  calculated.  Then  an  improved  estimate  of  SQ  is 


-  _■  R  -■  /II  \  R' 

S  =8„— - =  5  4-  —  - 

0  dR/dS'  0  V  6  M8(J  )2-l] 


(4^87) 


Table  3  shows  0g3  and  c2  tabulated  for  the  full  range  of  80.  Note  that  c2  changes 
only  slightly  and  therefore  would  not  be  an  accurate  means  of  reflecting  skewness. 
By  comparison,  the  relative  sensitivity  of  0g3  for  indicating  skewness  is  quite 
evident. 


TABLE  3  -  DEPENDENCE  OF  THIRD  MOMENT  AND 
SQUARE  OF  HALF  MAXIMUM  HEIGHT  ON  THE 
SKEW-BIAS  PARAMETER 


<$o 

9S3 

c2 

mm 

0 

9 

-0.009810 

9.002916 

wBm£. 

-0.019569 

9.011664 

0.054 

-0.029226 

9.026244 

0.072 

-0.038730 

9.046656 

0.090 

-0.048031 

9.072900 

0.108 

-0.057077 

9.104976 

0.126 

-0.065818 

9.142884 

0.144 

-0.074202 

9.186624 

0.162 

-0.082180 

9.246196 

0.180 

-0.089699 

9.291600 

0.198 

-0.096709 

9.352836 

0.216 

-0.103160 

9.419904 

0.235 

-0.108999 

9.492804 

0.252 

-0.114177 

9.571536 

0.270 

-0.118643 

9.656100 

0.288 

-0.122345 

9.746496 

0.306 

-0.125233 

9.842724 

0.324 

-0.127255 

9.944784 

0.342 

-0.128361 

10.052676 

0.3538132 

-0.128565 

10.125000 

CONCLUSIONS 


Knowledge  of  the  mixed  lubrication  process  and  the  associated  interactive 
wear  mechanisms  is  vital  in  the  engineering  of  a  variety  of  machine  components  for 
reliable  service  in  the  operational  environment  of  naval  ships.  The  thermo¬ 
mechanical  interactions  can  be  effectively  studied  by  computer-aided  analysis 
techniques.  However,  any  accurate  theoretical  model  must  penetrate  the 
traditional  boundaries  that  isolate  classical  disciplines  from  one  another.  The 
analytical  technique  presented  here  required  methodologies  to  be  combined  that 
were  developed  in  the  fields  of  numerical  computations,  statistics,  heat  transfer, 
continuum  mechanics,  and  topographic  analysis  of  surfaces,  in  addition  to  the 
theory  of  hydrodynamic  lubrication.  Experimental  verification  of  theoretical 
predictions  presented  will  be  necessary  to  prove  that  such  an  assemblage  of 
analytical  tools  is  valid,  so  that  these  tools  can  eventually  be  used  with 
confidence  in  the  engineering  of  actual  hardware. 


RECOMMENDATIONS 

Because  the  study  of  jnixed  lubrication  is  inherently  an  interdisciplinary 
undertaking,  a  well-focused  overview  must  be  balanced  by  commensurate  stimulation 
of  the  constituent  disciplines  which  play  the  pacing  roles.  At  present,  it  is 
timely  to  demonstrate  the  possibility  of  quantitative  analysis  of  the  mixed 
lubrication  process  in  the  realistic  operational  environment  of  a  machine 
component . 

A  recent  endeavor  (Appendix  A)  yielded  a  reasonable  analytical  prediction  of 
friction  behavior  of  the  tilt-pad  bearing  in  the  one-dimensional  approximation. 
Such  an  analysis  procedure  can  be  readily  upgraded  to  consider  a  full  complement 
tilt-pad  thrust  bearing.  Specifically,  we  recommend  that  a  computer-aided 
analysis  procedure  for  the  tilt-pad  thrust  bearing  be  developed  to  include  the 
following  features: 

o  variable  sliding  speed  due  to  radius  change, 
o  side  leakage  flows, 
o  surface  roughness  effects, 


o  film  cavitation  (Swif t-Stieber  condition), 

o  moment  balance  about  the  pivot  point,  and 

o  thermo-elastohydrodynamic  interactions  as  they  affect  temperature- 
dependent  fluid  viscosity  and  film  profile  distortions. 

We  concede  that  certain  aspects  of  surface  roughness  effects  can  be 
theoretically  modelled  only  rather  crudely  at  present.  Although  such  crude 
approximations  appear  adequate  to  describe  the  friction  behavior,  not  enough  is 
known  about  surface  deterioration  when  hydrodynamic  lubrication  is  incomplete. 
Traditional  wear  rate  estimation  based  on  empirical  test  data  is  suitable  only  for 
low  speed  machine  elements  for  which  the  relative  fit  and  conformance  to  a  nominal 
geometrical  profile  determine  the  functional  performance.  For  high  speed  machine 
elements,  lack  of  full  film  lubrication  usually  precipitates  violent 
thermomechanical  events  that  rapidly  lead  to  catastrophic  failure.  Our  goal  is  to 
control  processes  in  mixed  lubrication  so  that  inevitable  surface  deterioration 
proceeds  at  a  benign  pace.  A  quantitative  understanding  is  imperative  of  the 
closely  coupled  phenomena  of  surface  chemistry,  material  moduli!  and  hardness, 
local  heat  transfer,  material  behaviors  in  the  presence  of  severe  thermal 
gradients,  thermal  distortions  and  dynamics,  form  and  size  of  wear  debris,  and  the 
topographical  changes  occurring  at  the  sliding  contact,  for  these  phenomena 
interact  and  govern  the  capability  of  the  machine  element  to  endure  and  to  recover 
from  tribological  overload  that  exceeds  the  full  film  lubrication  capacity. 
Progress  in  this  area  can  be  fostered  by  a  well-coordinated  effort  which 
encompasses : 

o  friction  and  wear  tests, 

o  surface  mechanics,  including  superficial  and  near-surface  failure 
processes , 

o  topographical  analysis  for  corroboration  and  correlation  with  mixed- 
lubrication  analysis  and  wear  tests, 

o  testing  of  a  full  complement  thrust  bearing,  and 

o  the  formulation  of  a  cohesive  theoretical  model  of  the  evolutionary 
dynamics  of  topographical  changes. 


Because  of  the  possible  synergistic  interactions  of  various  processes,  the 
knowledge  of  mixed  lubrication  for  one  set  of  machine  components  and  environmental 
conditions  may  not  directly  apply  to  another  machine  component  or  to  substantially 
different  environmental  conditions,  even  though  the  same  fundamental  processes  are 
involved.  However,  the  general  approach  to  acquiring  the  knowledge  is  applicable 
to  many  other  machine  elements  that  are  vital  in  naval  ship  systems.  The  face 
seal  and  the  elaL  'ric  sterntube  bearings  are  candidate  topics  worthy  of 
attention  when  resources  are  available. 
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ABSTRACT 


The  term  "mixed  lubrication"  is  used  to  describe  regions  of  hydrodynamic 
lubrication  in  which  part  of  the  load  carried  by  a  bearing  results  from  actual 
metallic  contact  between  mating  bearing  surfaces.  Past  bearing  design  practice 
has  been  to  size  bearings  to  avoid  this  phenomenon.  However,  recent  interest  in 
reducing  size  and  weight  in  shipboard  machinery  has  led  to  a  reassessment  of  this 
mode  of  bearing  operation.  This  report  presents  a  computational  scheme,  based  on 
a  stochastic  model  of  roughness  distribution,  for  the  determination  of  bearing 
performance  in  the  mixed  regime.  The  model  also  considers  "crowning",  or  convex 
deformation  of  the  bearing  surface  that  arises  due  to  thermal  or  mechanical 
deformations  while  in  service.  This  is  known  to  cause  a  significant  rearrangement 
of  the  pressure  profile  on  the  bearing  surface.  Allowance  is  also  made  in  the 
determination  of  bearing  performance  for  the  frictional  moment  that  acts  on  a 
tilt-pad  bearing,  which  may  assume  sizeable  proportions  in  the  mixed  regime.  The 
results  are  presented  in  the  form  of  an  interactive  computer  program. 


I 


INTRODUCTION 


As  the  steam  engine  came  of  age  in  the  nineteenth  century,  the  need  for 
dependable  bearings  grew,  and  journal  bearings  were  successfully  developed.  But 
it  was  not  until  nearly  the  turn  of  the  century  that  the  scientific  theories  of 
Petrov,  Tower,  and  Reynolds  revealed  how  such  bearings  operated.  Prior  to  their 
efforts  it  was  believed  that  bearing  friction  depended  on  material  properties,  but 
theory  showed  that  this  was  not  the  case,  and  that  friction  was  merely  the 
consequence  of  viscous  shearing  stresses  developed  in  the  thin  lubricant  film  that 
separated  the  mating  bearing  surfaces. 

Application  of  these  theories  led  eventually  to  the  development  of  practical 
thrust  bearings.  Working  independently,  Kingsbury  in  the  United  States,  and 
Mitchell  in  Australia  invented  the  tilt-pad  thrust  bearing,  which  had  far  reaching 
consequences.  Among  other  things,  their  invention  greatly  facilitated  the 
conversion  of  the  world's  navies  and  commercial  fleets  from  sail  to  steam. 

In  the  intervening  years  a  vast  literature  dealing  with  various  aspects  of 
hydrodynamic  lubrication  has  emerged,  and  it  seems  safe  to  say  that  lubrication 
mechanics  is  well  understood,  at  least  when  viewed  from  the  standpoint  of  lightly 
loaded  bearings.  However,  in  recent  years  an  interest  has  developed  in  weight  and 
size  reduction  in  steam  propulsion  plants  designed  for  shipboard  use  dictating 
significant  increases  in  bearing  load  capacity  beyond  that  which  has  been 
previously  regarded  as  good  engineering  practice.  Such  bearings  typically  exhibit 
phenomena  that  cannot  be  explained  by  classical  hydrodynamic  lubrication  theory. 
Essentially,  at  higher  loads,  the  minimum  film  thickness  is  reduced  to  the  extent 
that  certain  amounts  of  metal-to-metal  contact  take  place  between  microscopic 
surface  irregularities,  or  asperities,  that  to  some  extent  exist  on  all  bearing 


surfaces.  This  is  known  as  "mixed  lubrication."  Interest  in  this  regime  stems 
from  the  fact  that  bearings  can  give  quite  satisfactory  performance  while 
operating  in  mixed  lubrication,  providing  that  the  load  carried  by  metallic 
contact  is  relatively  small  as  compared  with  the  total  load.  But  at  higher 
loadings,  complete  breakdown  of  the  film  occurs,  friction  rapidly  increases,  and 
the  bearing  fails.  Predicting  this  point  of  departure  is  a  primary  goal  of  mixed 
lubrication  theory.  Material  properties  play  an  important  role  in  this 
prediction,  indicating  that  earlier  theorizers  may  not  have  been  completely  wrong. 

While  the  breakdown  phenomenon  has  been  observed  experimentally  for  many 
years,  it  has  only  been  in  recent  times  that  serious  attempts  have  been  made  to 
model  this  process  mathematically.  Christensen!  has  presented  a  stochastic  theory 
of  mixed  lubrication  in  which  the  asperities  are  distributed  in  a  near-Gaussian 
way  and  examines  a  plane  slider  bearing  according  to  this  rationale.  His  results 
are  illuminating  in  that  the  model  gives  both  qualitative  and  quantitative 
performance  closely  matching  observed  behavior. 

Classical  bearing  theory,  when  applied  to  a  plain  slider  bearing,  also  fails 
to  account  for  the  operation  of  a  centrally  pivoted  bearing,  becoming 
mathematically  degenerate  for  this  case.  Despite  this  frustration,  centrally 
pivoted  bearings  give  fully  satisfactory  operation  and  are  in  widespread  use  in 
applications  requiring  thrust  to  be  carried  regardless  of  rotational  direction.  A 
good  explanation  of  this  evident  paradox  has  been  given  by  Abramovitz^,  who  showed 
that  an  exceedingly  small  amount  of  crowning,  or  convex  surface  curvature  on  a 
slider  face  (which  may  result  from  elastic  or  thermal  deformations,  or  polishing 
operations)  radically  alters  the  pressure  distribution  on  the  slider  in  a  way  that 
results  in  an  upstream  migration  of  the  pressure  center. 


COPY 


Finally,  we  call  attention  to  another  shortcoming  of  classical  theory  when 
applied  to  a  bearing  that  is  free  to  pivot,  and  that  is  the  essential  "static" 
nature  of  the  theory,  in  that  the  inclination  angle  of  the  bearing  pad  and  the 
ratio  of  leading  to  trailing  edge  film  heights  are  specified  a  priori  and  the 
bearing  is  then  analyzed  under  conditions  of  varying  load  and  speed.  It  is 
further  assumed  that  the  center  of  pressure  and  the  pivot  center  are  coincident, 
based  on  the  practical  assumption  that  the  offset  distance  of  the  center  of 
pressure  from  the  pivot  center  required  to  compensate  for  the  moment  due  to 
tangential  shearing  stresses  is  negligible.  While  this  simplification  appears  to 
be  fully  justified  in  the  pure  hydrodynamic  regime,  it  becomes  highly  suspect  in 
the  nixed  regime,  where  frictional  moments  can  reach  sizeable  proportions.  A  real 
bearing  pad  responds  to  variations  in  load  and  speed  within  its  operating 
capability  by  altering  its  inclination  angle  and  mean  film  height  so  that  the  load 
is  automatically  developed  and  the  moment  about  its  pivot  center  cancelled, 
regardless  of  operating  regime. 


1-1  Purpose  of  this  Study 

The  present  study  was  undertaken  to  explore  performance  boundaries  of  tilt 
pad  bearings  in  a  fashion  that  is  unencumbered  by  restrictions  of  classical  theory 
previously  discussed.  Considering  the  large  number  of  variables  in  this  class  of 
problem,  it  is  not  practicable  to  present  performance  in  customary  graphical 
(generalized)  formats.  In  the  final  analysis,  it  is  much  more  convenient  to 
compute  performance  for  each  specific  example  of  interest.  In  the  following,  we 
present  the  rationale  leading  to  a  computational  scheme  that  permits  numerical 


experimentation  for  a  class  of  curved  slider  bearings  operating  in  the  mixed 
lubrication  regime,  having  arbitrary  pivot  coordinates. 


II  APPROACH 


We  adopt  the  stochastic  model  of  Christensen  for  a  slider  bearing  free  to 
tilt  about  some  arbitrarily  specified  point.  The  macroscopic  film  height 
distribution  is  assumed  to  be  that  which  results  when  the  slider  surface  is 
parabolically  crowned,  having  a  maximum  crown  height  at  the  midpoint  of  the 
bearing  face.  Upon  the  macroscopic  height  distribution  is  superimposed  a  randomly 
distributed  series  of  linear  roughness  asperities,  aligned  so  that  the  peaks  and 
valleys  thus  formed  lie  in  the  direction  of  motion.  Thus,  at  any  point  along  the 
surface,  the  film  height  H  is  represented  as 

F  =  h(x)  +  h$  (x,z) 


(1) 


where  h  is  the  macroscopic  film  height  and  hs  is  a  randomly  distributed  roughness, 
or  asperity  height. 

In  regions  along  the  bearing  surface  where  the  film  height  is  greater  than 
the  asperity  height,  the  load  and  friction  are  assumed  to  be  purely  hydrodynamic 
in  nature,  whereas  in  regions  where  the  film  height  is  less  than  the  asperity 
height,  both  the  load  and  friction  are  carried  by  separate  contributions  of 
hydrodynamic  forces,  and  other  forces  that  depend  on  the  mechanical  properties  of 
the  bearing  material,  in  linear  combination.  Viscosity  is  assumed  constant,  and 
side  leakage  is  neglected.  Regarding  side  leakage,  it  should  be  noted  that  in 
regions  of  mixed  lubrication,  the  lubricant  is  effectively  trapped  in  the  valleys 
between  asperities,  a  fact  which  tends  to  mitigate  the  otherwise  sweeping 


assumption  of  no  side  leakage.  The  model  is  also  two-dimensional,  and  does  not 
include  centrifugal  effects. 

The  general  arrangement  and  nomenclature  are  indicated  in  Figure  1. 


Ill  FORMULATION 

III-l.  Basic  Relationships 

Under  the  stated  assumptions,  the  steady,  two-dimensional  form  of  the 
Reynolds  equation  is* 

(dJdx)  {E  (H3)  dp/dx)  =  6qt JdE(H)/dx 


(2) 


where  E(«)  is  the  expected  value  of  «,and  H  =  h(x)  +  hs(x,z),  the  hydrodynamic 
film  thickness.  The  value  of  p  in  the  above  expression  is  taken  to  be  the  mean 
hydrodynamic  pressure  at  any  location  x,  U  is  the  slider  speed  and  q  is  the 
viscosity  of  the  lubricant,  assumed  constant. 

The  random  roughness  distribution  is  assumed  to  be  well  represented  by  the 
near-Gaussian  function 

Khf)  =  (35132c?)  (c2  -  hs2)3 


where  2c  is  the  maximum  peak-to  valley  height  of  the  roughness.  Then  the 
expection  of  any  quantity  g(x)  is  obtained  from  the  expression 

E{g(x))  =  |  g{x)  f{h)dhg. 


(3) 


where  the  upper  limit  of  integration  is  always  c,  and  the  lower  limit  either  -c, 
or  -h,  depending  on  whether  the  film  thickness  is  greater  or  less  than  c.  This 
corresponds  to  regions  of  pure  hydrodynamic  and  mixed  lubrication.  In  particular, 
in  regions  of  pure  hydrodynamic  lubrication,  it  is  readily  shown  that 


E(K)  =  K,  a  constant 
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E(H)  =  h, 

and  E(H3)  =  h3  +  he2 13, 

whereas,  in  mixed  lubrication  regions,  more  detailed  calculations  show  that 

E(K)  =  (K/32c7)  {16c7  +  35c6 h  -  35c4 h3  +  21c2 h5  -  5h7}, 

E(H)  =  (35132c7  >(c8I8  +  (16l35)c7 h  +  (112k6 h2  -(H4)c4h4 
+  (U10)(?h6 -(V56)h8}, 


(4) 


and 

E(H3)  =  (35/32c7)  (1140k10  +  (16/I05)c9h  +  (3/8)c8h2 
+  (16l35)c7 h3  +  (l/4)c6h4  -  ( U20)c4h6 
+  (3l280)c?h6  -  (H840)h10}. 

The  nominal  film  height  h(x)  for  a  crowned  slider  bearing  (see  F-l)  is  given 

,  o 

by  the  expression 

h  =  h2(a  +  (l  -a)  (( l  +  4HJmB2)  x  -  (4/fc/mB3)x2)> 

(5) 


where  h2  is  the  trailing  edge  film  height,  a  =  h^/t^,  hj^  is  the  leading  edge  film 
height,  m  is  the  bearing  inclination  angle,  measured  in  radians  from  the  direction 
of  motion,  and  B  is  the  slider  length.  A  useful  relationship  between  a,  m,  B,  and 
h£  is 

h2  =  mB/  (a-l). 

(6) 


The  boundary  conditions  are  the  usual  specification  of  zero  gauge  pressures 
at  the  leading  and  trailing  edge  of  the  slider,  i.e, 

p(0)  =  p(B)  =  0. 
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111-2. 


Pressure  Distribution 


For  a  given  slider  bearing  characterized  by  its  physical  dimensions  (B,  Hc, 
c)  operating  at  a  speed  U  with  a  given  lubricant  viscosity  q,  it  is  possible  to 
obtain  a  solution  to  the  Reynolds  equation  by  direct  numerical  methods,  providing 
that  two  of  the  three  variables  a,  m,  and  h2  are  specified.  It  is  emphasized  that 
the  appropriate  expressions  for  the  expectancy  operators  must  be  employed  at  each 
point  in  the  integration,  depending  on  whether  the  lubrication  is  purely 
hydrodynamic,  or  mixed  at  each  point.  The  first  integration  is  straightforward, 
and  yields 


dpldx  =  6i\U 


(E(H)/E(H3)  +  C/E(H3) 


(7) 


where  C  is  a  constant  of  integration. 

The  second  integration  is  performed  by  numerically  calculating  the 
expressions 


p  =  Si\U 


X  r  x 

(. E(H)  IE{H3))  dx  +  C  (l  IE(H3))  dx 

o  '  o 


whereupon,  it  results  from  the  boundary  conditions  that 


■  B  B 

C=-  (EW/EiH3))  dxl  (1  /EiH^dx. 


III-3.  Bearing  Load  and  Center  of  Pressure 

Having  performed  the  above  integrations  which  give  the  pressure  distribution 
that  exists  on  the  slider  face,  taking  into  account  that  some  regions  are  in  mixed 
lubrication,  it  is  now  necessary  to  introduce  an  ad  hoc  assumption  regarding  the 
manner  in  which  the  total  load  is  developed.  The  simplest  model  of  the 
contribution  of  the  roughness  asperities  to  the  load  assumes  that  in  mixed 
lubrication  regions,  the  metallic  contact  load  is  equal  to  the  product  of  the 


yield  bearing  stress  of  the  weaker  bearing  material,  multiplied  by  the  fraction  of 
area,  per  unit  bearing  width,  over  which  such  contact  takes  place.  Thus,  the  load 
capacity  (per  unit  bearing  width)  is  represented  as 

2 


w 


fllA  f  *2  [B 

—  |  pdx  +  J  pdx  +  j  pdx, 


(s: 


where  in  the  first  and  last  integrals,  p  is  the  ordinary  hydrodynamic  pressure 
upstream  and  downstream  of  the  interference  region  that  exists  between  1^  and  l2> 
The  second  integral  is  the  sum  of  the  load  carried  by  ordinary  hydrodynamic 
pressure  in  the  valleys  between  asperity  peaks,  and  the  load  carried  by  the 
asperities;  i.e., 

[l2  f  *2  f  *2 

I  pdx  =  I  r(x)pdx+py  I  (1  —  r)dx 

1  ll  J,1 

o: 


where  Py  is  the  yield  stress  of  the  weaker  bearing  material,  and  r  is  the  fraction 
of  bearing  area  formed  by  the  valleys  between  asperity  elements,  calculated  from 
equation  (4)  to  be 

r  =  (1132c7)  {16c7  +  35c6 h  -  35c4 h3  +  21c2 h5  -  5h7}, 


per  unit  bearing  width. 

The  center  of  pressure  is  obtained  in  a  straightforward  fashion  by  weighting  the 
above  integrands  with  x,  performing  the  indicated  integrations,  and  dividing  by 
the  load;  i.e.. 


V  l„  pxdx/l0  pdx • 

(10 


where  it  is  to  be  '«r  ood  that  the  integrations  follow  the  above  rationale 


1 1 1-4 


Slider  Friction 


The  shearing  stress  acting  on  each  element  of  the  slider  surface  is  obtained 
from  the  expression* 

t(x)  =  qUE(l/H)  +  1/2  (dpldx)E(H), 

(11 

where  dp/dx  is  given  by  (7). 

The  total  frictional  force  on  the  slider  is  determined  in  a  manner  similar  to 
the  method  used  in  obtaining  the  load;  i.e., 

*2 


f1!  f  *2  [B 

=  idx  +  idx  +  tdx, 
Jo  J 1.  Jl„ 


(12 


where,  here  again. 


f  ^  =('2 


n  dx  4-  o  (1  —  r)dx, 
1  Jl. 


(13 


and  o  is  the  yield  shearing  stress  of  the  weaker  bearing  material. 

The  frictional  moment  about  the  pivot  center  is  the  integral  of  the  shear 
acting  on  each  surface  element  weighted  by  its  distance  from  the  y  coordinate  of 
the  pivot  center.  Since  the  bearing  surface  is  curved,  this  distance  is  not, 
strictly  speaking,  constant  along  the  bearing  surface.  However,  in  our  area  of 
interest,  the  maximum  ratio  of  crown  height  to  slider  length  is  of  the  order  10~3 
and  it  is  sufficiently  accurate  to  assume  that  all  shearing  stresses  act  at 
constant  distance  from  the  y  pivot  coordinate,  thus  simplifying  subsequent  moment 


balance  calculations 


In  order  to  perform  the  indicated  integrations  to  obtain  the  frictional 
force,  the  additional  expectation  function  E(l/H)  must  be  evaluated.  In  the  pure 
hydrodynamic  region,  this  is,  from  (3), 

£(1/H)  =  (35/32c7)  (c2  -  A2)3  1  n{(H  +  cV(A-c)} 

+  (2/ 1 5  )cA(  1 5  /i4  ~40c2A2  +  33c4) 

and,  in  the  mixed  region,  the  equivalent  expression  is 

E(MH)  =  35/32c7)  (c2  -  h2)3  In  {(A  +  c)/(h  -  A*)}  +  (8A3/3)(A  +  c)3 

-(1/6)  ( A  +  c)6  —  4(c2  -  h2)h{h  +  c)3 
4-  6(c2  -  h2)2h(h  +  c)  +  6(c2  -  h2)h2)(c  +  h2) 

—  (1/6)(A  +  c)6  —  4(c2  -  h2)h{h  +  c)3 
—  3A2(c  +  hf  +  (6/5 )h(h  +  cf  +  (3/4) (c2  -  h2)(h  +  c)4 

—  (3/2)(c2  —  A2)2  (A  +  c)2 

The  term  h*  in  the  above  expression  appears  in  order  to  prevent  the 
logarithmic  term  from  becoming  unbounded.  Following  Christensen,  we  set  h  -  h*  = 
mB/1000  in  subsequent  calculations.  Numerical  experimentation  shows  that  the 
resulting  solutions  are  not  particularly  Sensitive  to  this  assumption. 
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III-5«  Coefficient  of  Friction 

The  frictional  force  indicated  above  acts  on  the  slider  surface.  The  total 
resistance  to  motion  is  this  force,  plus  the  component  of  load  acting  in  the  plane 
of  rotation*  A  precise  determination  of  this  second  contribution  would  involve  an 
incef'  *ation  of  the  components  of  load  in  the  direction  of  rotation  that  act  on 
each  element  of  the  curved  slider  surface.  Since,  however,  the  radius  of 
curvature  of  the  curved  surface  is  many  orders  of  magnitude  greater  than  the 
slider  length,  it  is  sufficiently  accurate  to  write 

F  =  F  +  mw, 

r  $ 

where  Fr  is  the  frictional  force  on  the  runner.  Then  the  coefficient  of  friction 
p  is  obtained  by  dividing  Fr  by  the  load;  i.e., 

p  =  FJw  +  m. 

I I 1-6.  Moment  Balance 

Referring  to  Figure  2,  it  is  seen  by  inspection  that  the  condition  that 
satisfies  zero  unbalanced  moments  about  the  pivot  center  is 

FsYP  =  <Xp-Cp)^,or 
Cp  =  *p-(FJw)Yp 

(14 

where  Xp,  Yp  are  fixed  by  design,  and  Cp,  Fg,  and  w  are  functions  of  speed,  load, 
etc.  the  term  Fg/w  is  recognized  from  the  above  as  a  first  approximation  to  the 
friction  coefficient  which,  for  bearings  operating  in  the  pure  hydrodynamic 

•  .  _  'l  ,  , 

regimes,  is  a  number  approximating  10  .  In  classical  formulations,  the  second 

term  is  neglected.  However,  in  the  mixed  regime,  this  simplification  is  not 
generally  justified,  and  thus  will  not  be  adopted  in  this  study. 

This  completes  the  basic  formulation  of  the  problem. 


IV 


COMPUTATION  OF  SPECIFIC  BEARING  PERFORMANCE 


IV-1«  The  Inverse  Problem 

Discounting  obvious  embellishments  required  to  accommodate  mixed  lubrication 
and  crowning  considerations,  the  preceding  analysis  has  followed  the  classical 
approach  to  bearing  performance;  i.e.,  the  Reynolds  equation  is  integrated  to 
obtain  the  pressure  distribution,  which  is  then  successively  integrated  to  obtain 
the  load  and  pressure  center.  A  parallel  integration  gives  the  frictional  force. 
Application  of  this  approach  depends  upon  practical  assumptions  regarding  the  mean 
film  height  (or  the  film  height  ratio,  which  amounts  to  the  same  thing)  and  the 
slider  inclination  angle.  When  this  is  done,  the  performance  of  a  given  bearing 
may  be  determined.  In  particular,  the  load  and  pressure  center  of  a  bearing 
operating  under  given  conditions  may  be  obtained  from  explicit  relations  (or,  in 
the  more  complex  case  treated  here,  by  direct  numerical  integrations)  that  take 
the  form 

w  =  w(a,m,U,q,  dimensions,  material  properties) 

Cp  =  Cp  (a,  m,U,q,  dimensions,  material  properties) 
where  a  and  m  are  assumed  known  to  reasonable  approximation. 

In  the  real  situation,  however,  the  pivot  point  is  fixed  by  design,  and  the 
load  and  speed  are  independent  variables.  As  these  vary  (often  simultaneously, 
according  to  some  speed-load  schedule)  the  film  height  and  inclination  angle 
automatically  adjust  to  accommodate  the  new  operating  conditions;  we  may  thus 
think  of  a  bearing  as  a  computer  that  continuously  solves  the  inverse  problem 
according  to  functional  relations  of  the  form 

a  =  a(w,  X '  ,  U,  etc.) 
m  =  m(u>,  Xp,  etc.) 
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where  Xp  and  w  are  specified. 

The  method  by  which  the  bearing  performs  this  inversion  cannot  be  duplicated 
by  rational  analysis.  Thus,  if  a  bearing's  performance  is  to  be  determined  across 
a  spectrum  of  operating  conditions,  it  must  be  done  by  trial-and-error ,  using  the 
forward  integration  technique  indicated  above.  It  is  evident  that  a  dual 
iteration  must  be  performed  on  the  variables  a  and  m,  which,  when  acting  in  proper 
(and,  one  would  hope,  unique)  combination,  results  in  a  correct  load  prediction, 
and  at  the  same  time  gives  assurance  that  the  dual  moments  about  the  pivot  center 
caused  by  load  and  friction  will  be  effectively  canceled. 

Computational  methods  exist,  of  course,  for  dual  iteration  problems. 

However,  such  schemes  usually  owe  their  success  to  the  fact  that  the  values  sought 
are  single  valued,  and  vary  monotonically  with  the  iteration  variables. 

Considering  the  sizeable  number  of  variables  in  the  present  problem,  there  is  no 
assurance  whatever  that  these  conditions  can  always  be  met.  It  is  clear  that 
physical  insights,  leading  to  "intelligent  guesses,"  must  be  sought. 

IV-2.  Inverting  the  Plane  Slider  Bearing  Problem 

By  way  of  illustrating  the  above  considerations,  as  well  as  obtaining  a 
highly  useful  means  of  arriving  at  an  intelligent  first  guess  needed  to  expedite 
the  iteration  process  in  the  more  complex  case  treated  in  this  study,  the  plane 
slider  bearing  problem  is  inverted;  i.e.,  given  the  load,  speed,  dimensions, 
viscosity,  and  pivot  center,  the  problem  is  to  obtain  the  inclination  angle  and 
film  height  ratio,  for  the  case  of  the  smooth  plane  slider  bearing. 

To  begin,  the  Reynolds  equation  is  solved  for  the  pressure  distribution.  The 
result,  expressed  in  partial  fractions  suitable  for  further  integration,  is: 
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(a-l)A2  1  l+ta-lk/B  (a  +  l)(l+(a-l)  (x/B)2  <a  +  l>  I 
which,  when  integrated  on  x  between  the  limits  0,  B,  gives 

6 nUB2  f,  2(a—  1)  I 

w  =  - —  j  ln(a)  -  - > 

(a-l)2/i2l  (a  +  1 )  I 

With  the  help  of  (6),  this  may  be  rewritten 

6 nU  [,  ,  t  2(a  —  1) 

w  =  — —  ln(a)  -  - 

m2  l  (a  +  1) 


=  6fi(/g(a)/m2 


The  center  of  pressure  is  obtained  by  a  second  integration  to  yield 

2a(a  +  2 )ln(a)  -  (5a  + 1 )  (a  - 1 ) 


2(a  — l)[(a  +  l)m(a)  -  2(a-l)| 
=  cp(a). 


Now  clearly,  if  it  were  possible  to  invert  this  last  expression,  yielding 

a  =  a(cp) 

then  m  would  immediately  follow  from  (16),  i.e., 

m  =  y/6r\U  f[a)!w. 


Unfortunately,  the  transcendental  form  of  (17)  prohibits  a  direct  inversion; 
however,  for  reasons  that  are  not  readily  apparent,  an  excellent  correlation 
formula,  good  to  1Z  of  the  exact  prediction  given  by  (17)  in  the  range  1  <.  a  .<8 


cp  =  .0919  ln(a)  +  .504 


which  is  readily  inverted  to  yield 


a  -  exp 
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10.88(cp-  .504)  | 

(2 

Having  now  determined  a,  the  value  of  m  follows  from  (18),  and  h2  from  (6).  The 

slider  friction  per  unit  bearing  width  is  obtained  from  the  expression 

_  nUB  6  2  ln(a) 

s  h2  a  4-  1  (a  - 1) 

the  runner  friction  is 

Fr  =  F3  +  mw, 

and  the  coefficient  of  friction  is  Fr/w.  Thus,  the  problem  is  inverted.  Note 
that  if  the  pivot  coordinates  are  specified,  it  is  possible  to  compute  the  actual 
pressure  center  required  to  nullify  the  load  and  frictional  moments.  Then,  with 
the  aid  of  (20),  a  new  value  of  a  is  found,  and  the  calculations  are  repeated. 
Study  shows  that  one,  or  at  the  most  two,  iterations  is  all  that  is  required  to 
completely  specify  the  frictional  force,  with  full  assurance  that  moment  balance 
is  achieved. 

There  is  an  interesting  extension  to  this  development  that  bears  directly  on 
the  more  complex  cases  being  dealt  with  in  this  study.  Returning  to  the 
expression  for  the  slider  friction  given  above,  then,  with  the  help  of  (6),  this 
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Thus  the  coefficient  of  friction  is,  again  using  (16), 


p  =  m 


f  (a)  +  6g(q) 


/6  g(a)  =  mF(a). 


(21 


The  full  solution  for  this  problem,  given  the  load,  speed,  viscosity,  and 
pivot  center  is  thus 

_  2(q+  l)/n(a)  -  3(q-l) 

11  ~  3(q  +  1  )ln(a)  -  6(q-l) 


w 


ln(a )  - 


2(q-l) 


(q  +  1) 


(22 


where  a  is  obtained  from  the  correlation  formula 


q  =  exp 


10.88(c  -  .504) 

p 


and  the  minimum  film  height  is 


_  mB 
2~  (a  —  1 ) 


w(a  - 1 ) 


ln(a)  — 


=T?II 


Note  that  (22)  may  be  written 


'‘=Cfa)(f )  • 


(20 


(6,  23 


(24! 


expressing  the  well-known  fact  that  the  coefficient  of  friction  is  proportional  to 
the  product  of  the  speed  and  viscosity,  divided  by  the  load,  all  raised  to  the 
one-half  power.  Once  again,  if  moment  balance  is  considered  important,  the  above 
relations  offer  a  rapid  iteration  method. 


IV-3  Effects  of  Mixed  Lubrication 

In  the  preceding  development,  a  simple  strategy  emerged  for  determining  the 


friction  coefficient  while  requiring  moment  balance.  In  the  mixed  lubrication 


regime,  there  is  no  hope  of  duplicating  this  since  more  variables  are  introduced, 
and  the  equations  do  not  admit  closed  form  solutions.  Thus,  the  problem  of 
finding  both  the  inclination  angle  and  the  ratio  of  film  heights  that  deliver  a 
given  load,  while  assuring  moment  balance,  is  formidable. 

But  if  we  pause  to  consider  that  in  general  the  primary  interest  is  not  so 
much  in  the  determination  of  these  parameters  for  a  given  load,  but  rather  a  range 
of  loads,  then  the  picture  is  somewhat  brighter.  In  the  pure  hydrodynamic  regime, 
(6)  shows  that  for  a  given  bearing  operating  at  a  given  speed,  the  inclination 
angle  and  load  exist  in  simple  relationship.  Thus,  given  the  load,  m  is 
established,  but  one  could  just  as  well  specify  m  and  compute  the  load  that  this 
corresponds  to.  Therefore,  in  a  numerical  determination  of  bearing  performance,  a 
reasonable  approach  seems  to  be  one  in  which  the  inclination  angle  is  successively 
stepped,  with  film  height  iterations  that  produce  moment  balance  being  performed 
at  every  step  along  the  way.  This  is  in  fact  the  approach  eventually  adopted 
after  considerable  experimentation  with  dual  iteration  schemes. 


V  NONDIMENSIONAL  SYSTEM 


In  order  to  facilitate  computations  and  to  illuminate  certain  consequences  of 
mixed  lubrication,  it  is  convenient  to  normalize  the  relevant  equations  of  motion. 
To  do  this,  we  note  that 

U  =  LIT,  p  =  FIL2,  x  =  L,  E(Hn)  =  Ln, 
where  F,  L  indicate  the  dimensions  of  force  and  length. 

Defining 


x '  =  x/B,  H’  =  H/h2,  p  '  =  php  qUB,  then 
d/dx  =  d/Bdx',  and  the  Reynolds  equation  (2)  becomes 

dE(H') 


—  EM'3) 
dx' 


=  6 


d_ 

dx' 


103 


dx' 


I  rJMMMBMWa 


COPY 


which  is  now  fully  nondimens ional . 

The  macroscopic  height  expression  (5)  may  be  normalized  on  h2  and  B  to  yield 


h'  =  hJh2=a-  (1 


-a)]  (1 


+  4H'/m)x'-  (4H'/m)x'£  , 

C  C 


where  H,e  =  Hc/B. 


Likewise,  all  expectation  functions  can  be  normalized  on  h2,  replacing  c  with 
c'  =  c/h2andh'  =  h/h2-  Then  it  follows  that 

E(Hn)  =  h2n  E(H'n),  where  n  is  any  integer. 

Following  this  rationale,  the  load,  given  by  (8)  may  be  rewritten 

r  f.  r  r2  r  >2  r  I* 

w  =  nUB2/hi  p'dx'+  I  rp'dx' +  (p  hi  /  nUB2)  \  (l-r)dx'  +  I  p'dx’ 

U°  Jvt  >2  Jr,  Jv2 

since  r'  =  r  from  above. 

The  nondimensional  shear  stress  i'  =  th2/nU  comes  from  (11)  and  is  given  by 

r’  =  E(l/H')+  midp’IdxWH') 
and  the  force  on  the  slider  is  therefore  written 

Fg=  nUB/h^  J  1  i'dx'  +  J  2  i’rdx'  +  (oh^qUB)  |  2{l-r)dx'  +  j  t'dx' 

With  this  nondimensionalization,  it  is  possible  to  generalize  and  thereby 
greatly  facilitate  numerical  integration  of  the  above  equation. 

Another  result  of  this  development  is  the  indication  that  bearing  performance 
cannot  be  predicted  by  the  simple  relationship  given  by  (24)  in  the  mixed  regime, 
owning  to  the  introduction  of  the  new  parameters  Py  and  o. 


VI  COMPUTATIONAL  SCHEME 

On  the  basis  of  the  foregoing,  a  computer  program  was  developed  in  the  BASIC 
language  for  use  on  an  interactive  system  such  as  the  Tektronics  4052.  The  inputs 
required  are  the  bearing  dimensions,  rotational  speed,  asperity  and  crown  heights, 
the  bearing  material  properties,  and  the  pivot  coordinates,  all  in  the  dimensions 


>  -f  .»  ^  - 
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indicated  in  the  program.  The  program  then  calculates  a  reasonable  first 
approximation  to  the  value  of  a  (=  hi/h2)  and  a  slope  m,  based  on  rationale 
previously  presented.  The  equations  are  then  successively  integrated  to  obtain 
the  friction  and  pressure  moments.  If  the  net  unbalanced  moment  is  greater  than 
the  load  multiplied  by  one-tenth  of  one  percent  of  the  bearing  length,  a  new  value 
of  a  is  tried,  calculated  from  a  "hunting"  routine,  and  the  process  repeated  until 
the  moment  balance  criterion  is  met.  Regarding  the  adoption  of  this  criterion,  it 
is  noted  that  this  particular  standard  appeared  to  represent  the  best  compromise 
between  reasonable  moment  balance  and  iteration  time,  which  can  be  quite  lengthy, 
especially  in  the  mixed  regime. 

When  the  moment  balance  criterion  has  been  met,  the  operator  then  receives  a 
listing  of  pertinent  bearing  parameters  such  as  the  load  and  friction  forces 
(broken  down  into  contributions  from  pure  hydrodynamic,  mixed  hydrodynamic,  and 
asperity  contact  forces),  the  shear  and  pressure  moments,  and  so  on.  The  operator 
may  then  select  a  new  slope,  and  the  process  is  then  repeated,.  Operators  who 
wish  to  change  the  initial  inputs  may  access  this  routine  by  inputting  a  zero 
slope.  Instructions  are  given  throughout  the  program  to  assist  the  operator  in 
selecting  workable  combinations.  For  example,  if  zero  crown  height  is  selected, 
and  the  pivot  center  is  specified  to  lie  at  the  midpoint  of  the  bearing,  the 
operator  will  be  informed  that  this  is  an  unworkable  combination  of  parameters, 
and  will  be  given  the  opportunity  to  try  new,  workable  combinations. 

As  shown  in  Ref.  2,  and  borne  out  in  this  study,  certain  regions  of  crowned 


bearing  operations  lead  to  negative  (gauge)  pressures  at  some  points  on  the  slider 
face,  which  may  in  some  applications  result  in  cavitation,  thus  altering  the 
pressure  distribution.  It  was  not,  however,  possible  to  pursue  this  aspect  within 


the  scope  of  this  study,  and  so  it  was  decided  to  merely  alert  the  operator  to 
this  possibility  when  negative  pressures  appear  through  the  statement,  "Bearing 
may  cavitate."  When  this  statement  appears,  the  operator  should  exercise  some 
caution  in  applying  the  results. 

The  requirement  for  continuous  slope  information  requires  continuous  inputs 
from  the  operator.  In  an  earlier  version  of  this  program,  the  slope  was  stepped 
successively  through  the  program,  but  it  was  found  that  very  small  changes  in 
slope  bring  about  quite  large  variations  in  load  and  friction  coefficient  in  the 
transition  regime  between  pure  hydrodynamic  and  mixed  lubrication,  and  this  scheme 
had  to  be  abandoned  because  too  much  information  was  being  "lost"  with  this 
routine. 

As  a  matter  of  information  to  operators,  it  is  noted  that  the  program 
"remembers"  the  value  of  a  that  was  found  by  iteration  in  a  previous  run,  and  uses 
this  as  the  first  estimate  for  the  next  run.  Therefore,  iterations  are  minimized 
if  the  operator  takes  relatively  small  "bites"  in  specifying  new  slopes.  That  is, 
the  larger  the  increment  in  m  for  each  succeeding  run,  the  greater  the  number  of 
iterations  will  be  required. 

A  complete  program  Listing  is  given  in  Appendix  I. 

VII  APPLICATION  OP  RESULTS 

To  illustrate  the  use  of  this  effort,  a  bearing  pad  presently  under  study  by 
the  David  W.  Taylor  Naval  Research  and  Development  Center  was  used  as  a  prototype. 
The  general  characteristics  of  this  slider  are: 

B  =  5. 87cm  (2.31  in.),  mean  diameter  =  13.2cm  (5.2  in.),  q  =  4.13  poise  (6  x  10-6 
reyns),  Py  =  1.22  x  10®  N/m^  (18,000  psi),  o/Py  =  0.2,  Xp  =  0.5,  yp  =  0.3. 
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With  these  fixed  inputs,  several  cases  of  asperity  height,  crown  heights,  and 


rotational  speeds  were  chosen  for  study.  These  are  summarized  in  the  following 
table: 


Case 

Hc 

c 

RPM 

I 

0.0002  in. 

0.0002  in. 

35 

II 

0.0002  in. 

0.0001  in. 

35 

III 

0.0002  in. 

0.0002  in. 

20 

IV 

0.0002  in. 

0.0002  in. 

50 

Figure  3  shows  the  effect  of  asperity  height.  Characteristic  of  this  class 
of  problems,  extremely  small  dimensional  changes,  as  compared  with  the  macroscopic 
dimensions  of  the  bearing  pad,  lead  to  significant  differences  in  the  performance 
characteristic.  This  is  reinforced  in  Figure  4,  which  shows  the  ratio  of  film 
heights  plotted  against  the  speed-load  parameter. 

Figure  5  shows  the  "nonsimilarity”  of  solutions  in  the  mixed  lubrication 
regime,  as  discussed  in  Section  V.  Finally,  Figure  6  shows  the  level  of  agreement 
between  the  calculated  values  of  shear  and  pressure  moments,  as  discussed  in 
Section  VI.  Here,  it  can  be  seen  that  the  level  of  agreement  increases  at  the 
higher  loads,  coresponding  to  the  moments  developed  in  the  mixed  regime,  where 
such  agreement  is  most  important. 

Many  similar  cases  were  explored  during  the  course  of  this  study,  leading  to 
the  same  essential  conclusions.  It  is  emphasized  that  the  utility  of  this  study 
cannot  be  realized  without  the  aid  of  a  computer,  which  allows  a  plethora  of  such 
examples  to  be  considered. 


VIII  SUMMARY 


This  study  represents  an  extension  of  prior  work  by  adapting  the  stochastic 
model  of  mixed  lubrication  due  to  Christensen  to  the  problem  of  the  curved  slider 
bearing  treated  by  Abramovitz.  As  a  further  refinement,  the  bearing  pad  is  left 
free  to  seek  the  combination  of  parameters  that  causes  the  dual  moments  due  to 
load  and  friction  to  balance,  within  a  narrow  margin.  The  results  are  presented 
in  the  form  of  an  interactive  computer  program  which  offers  the  operator  a  wide 
range  of  parametric  authority  in  the  design  or  analysis  of  a  family  of  roughened, 
crowned,  tilt  pad  bearings. 

As  Christensen  points  out,  no  allowance  is  made  for  wear  of  the  asperity 
contacts  in  the  mixed  regime,  which  almost  certainly  would  exist  to  some  greater 
or  lesser  extent  in  the  real  situation.  The  major  consequence  of  this  restriction 
is  probably  an  overly-pessimistic  prediction  of  the  friction  coefficient  at  the 
higher  loads.  Thus,  the  error  introduced  by  this  assumption  is  on  the  side  of 
conservatism. 

It  may  also  be  noted  that  the  crown  height  and  viscosity  are  considered 
constant  in  this  model  (although  they  may  easily  be  changed  from  one  run  to 
another)  whereas  in  the  real  situation  they  are  also  most  probably  time-dependent 
due  to  thermal  effects.  This  effect  cannot  be  accounted  for,  of  course,  unless 
heat  transfer  is  introduced  into  the  present  model.  This  would  seem  to  be  the 
logical  next  extension  of  the  present  work. 
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100  PRINT  "DO  YOU  WANT  AN  EXPLANATION  OF  THIS  PROGRAM?” 

110  INPUT  P$ 

120  IFSEG$(P$,  1, 1)  =  "Y"THEN  3120 

210  PRINT  "INPUT  RPM” 

220  INPUT  U1 

230  PRINT  "INPUT  MEAN  BEARING  SURFACE  DIAMETER,  INCHES” 

240  INPUT  K1 

250  PRINT  "INPUT  VISCOSITY,  REYNS  X  1E6” 

260  INPUT  M2 
270  M2  =  M2*lE-6 

280  PRINT  "INPUT  SHOE  LENGTH,  INCHES” 

290  INPUT  B 

300  PRINT  "INPUT  ASPERITY  HEIGHT,  INCHES” 

310  INPUT  C2 

320  PRINT  "INPUT  YIELD  BEARING  STRESS,  PSI  X  IE-3” 

330  INPUT  P8 

340  P8=P8*1E3 

350  PRINT  "INPUT  RATIO  OF  SHEAR  TO  YIELD  STRESS” 

360  INPUT  M3 

370  PRINT  "INPUT  CROWN  HEIGHT,  INCHES” 

380  INPUT  H4 

390  PRINT  "INPUT  PIVOT  CENTER,  X,  Y,  FRACTION  OF  SHOE  LENGTH” 

400  INPUT  P2,  P3 

410  U2=U1*(3. 1416/60)*K1 

440 

450  K5  =  0 

460  A=  0042*EXP(10.88*P2) 

470  A  =  A  +  (H4/B>*1E4 

480  IFA>1  THEN  500 
490  A=  1.001 
500  Z1  = .  005 
510  M5  =  .00025*(P2  +  .5) 

520  K9  =  0 

530  K8=0 

540  N8=0 

560  K7  =  0 

570  PRINT 
580  LI  =  0 
590  16=0 

600  17  =  0 

605  X4= 0 

606  X5  =  0 

610  18=0 
620  J1  =  0 
630  J6  =  0 
640  J7  =  0 
650  J8  =  0 
660  N9  =  0 
670  11  =  0 

680  12  =  0 
690  13  =  0 
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14=0 
X4=0 
X4=0 
T1  =  0 
T2  =  0 
T3  =  0 

H2  =  M5*B/(A-1) 

M6=H2/B 
H7  =  H4/B 
X5  =  0 

G=P8*M6A2/(M2*U2) 

H  =  G*M3/M6 
Ml  =  M5*B/1000 
Ml  =M1/H2 

FORX=  Z1 T01  STEP  Z1 
F  =  H4/H2 

FI  =  A  +  (1-A)*((1  +  4*H4/(M5*B))*X-(4*H4/(M5*B))*X  *  2) 

C1  =  C2/H2 

IF  Cl/Fl  <  1  THE  N  1000 

Nl  =  (Cl*8)/8  +  (16/35)*C1*(C1*6)*F1  +  5*F1*(C1*6)*F1 
N1  =  N1  -(.25)*(C1*F1)*4  +  (C1*2)*(F1*6)*  l-(Fl*8)/56 
N1  =  N  1*35/(32*C1 '7) 

D2  =  (Cl‘10)/40  +  (16/105)*F1*C1*9  +  (3/8)*(Cl  *8)*(F1‘2) 

D2  =  D2  +  (16/35)*(C1*7)*F1*3 

D2  =  D2  +  25*(C  1*  6)*F1  *  4  -  ( 1/20)*C  1*  4*F1 *  6  +  3/280 

D2  =  D2-3/280  +  <3/280)*(Cl  *2)*(F1  *8)-(  1/840)*(F1  *  10) 

D2  =  35*  D2/(32*C  1*7) 

D1  =  D2/N1 
GOTO  1020 

D1  =  F1  *  (FI  f  2  +  Cl  T  2/3) 

D2=  F1*D1 

11  =  11+  Z1/(D1) 

12  =  I2  +  Z1/(D2) 

NEXT  X 

C  =  11/12 

FOR  X  =Z1  TO  1  STEP  Z1 

FI  =  A  +  (1-A)*((1  +  4*H4/(M5*B))*X-(4*H4/(M5*B))*X  *2) 

Cl  =  C2/H2 

IF  Cl/Fl  C1THEN  1460 
L1  =  L1  +Z1 
N9  =  N9  + 1 
X4=X  +  Z1*(1-N9) 

R1  =  (1/(32*C1  *  7))*(16*C1  *  7  +  35*(C1  *  6)*F1  -35*(C1  *  4)*F1  *  3) 
R1  =  R1  +  (1/(32*C1  *7))*(21*(C1  *2)*F1  *5  -5*F1  *7) 

I6  =  I6  +  (1-R1)*Z1 
J6=J6  +  (1-R1)*X*Z1 

N1  =  (Cl  *8)/8+(16/35)*(Cl  *7)*F1  +  5*(C1  ‘6)*(F1  *2) 

N1  =  Nl-  .25*(C1‘4)*F1‘4  +  .1*(C1‘2)*(F1  *6M1/86)*(F1)’8 
N1  =  35*N1/(32*C1*7) 

D2  =  (Cri0)/40  +  (16/105)*(Cl’9)*Fl  +(3/8)*(Cl‘8)*Fl*2 
D2  =  D2+  ( 16/35)*(C  1 "  7)*(F1  *  3)  + .  25*(C1  *6)*F1  4 
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1220  02  =  D2-.05*(C  1‘4)*(F1*6)  +  (3/280)*(C  1  *  2)*(F1  * 8)-( 1/840)*(F1  *  1 0) 

1230  D2  =  35*  D2/(32*C  1*7) 

1240  D1  =  D2/N1 

1250  N3  =  ((Cl  “2-F1  *2)*3)*L0G((F1  +  Cl)/Ml)  +  (8*F1‘3/3)*(F1  +  Cl)‘3 

1260  N3=N3-(1/6)*F1  +  C1)*6-4*(C1*2-F1‘2)*F1*(C1  +  Fl)*3 

1270  N3  =  N3  +  6*(C1‘ 2-F1  ‘2)'2*F1*(F1  +  Cl)  +  6*(C1  2-F1  ‘2)*F1  *2*(Cl  +  FI)*  2 

1280  N3  =  N3-3*(Fl*2)*(Cl  +  Fl)‘4+(6/5)*Fl*(Cl  +  Fl)*5 

1290  N3  =  N3  +  (3/4)*(Cl  *2-Fl  *2)*((F1  +  Cl)*4) 

1300  N3  =  N3-(3/2)*(Cl  *2-Fl  *2)*(Cl  *2-Fl  *2)*(F1  +  Cl)‘2 
1310  N3  =  N3*(35/(32*C1  *7)) 

1320  T1 =T1  +  R1*Z1*(N3  +  3*N1*(N1-C)/D2) 

1330  I3  =  I3  +  Z1/D1 

1340  I4=I4  +  Z1/D2 

1350  Q  =  I3-C*I4 

1360  IF  Q>-lE-4  THEN  1410 

1370  N8=N8  +  1 

1410  Q=Q*6 

1420  I7  =  I7  +  Q*Z1*R1 

1430  J7  =  Q*Z1*R1*X  +  J7 

1440  X5  =  X 

1450  GOTO  1710 

1460  D1  =  FI  *  (Fl  t  2  +  Cl  f  2/3) 

1470  D2=D1*F1 

1480  Nl  =  Fl 

1520  IF  Cl  <  IE-3  THEN  1570 

1521  IF  F1>0  THEN  1530 

1522  PRINT  "YOU  WILL  HAVE  TO  INCREASE  SLOPE  FOR  A  SOLUTION” 

1523  GO  TO  2506 

1530  N3  =  ((Cl'2-Fl‘2)‘3)*LOG((Fl  +  Cl)/(F1-Cl)) 

1540  N3=N3+  (2/15)*(Cl*Fl)*(15*Fl"4-40*(Cl*Fl)*2  +  33*Cl‘4) 

1550  N3  =  N3*(35/(37*C1‘7)) 

1560  GOTO  1580 
1570  N3=  1/FI 

1580  I3  =  I3  +  Z1/(D1) 

1590  14  =  14  +Z1/(D2) 

1600  Q=I3-C*I4 

1610  IF  Q>-lE-4  THEN  1660 

1620  N8=  1 

1660  LETQ  =  Q*6 

1670  N3=  1/FI 

1680  T2  =  T2  +  (3*N1*(N1-C)/D2  +  N3)*Z1 

1690  I8  =  I8+Z1*Q 

1700  J8=J8  +  Z1*Q*X 

1710  NEXT  X 

1720  W  =  I7  +  I8  +  G*I6 

1730  W1  =  W*U2*M2/M6'2 

1740  W2  =  I7*U2*M2/M6  2 

1750  W3  =  I8*U2*M2/M6‘2 

1760  W4=G*I6*U2*M2/M6  2 

1780  T3  =  Tl  +  T2  +  I6*H 

1790  T5  =  T2*U  2*  M2/M6 
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1800  T6  =  T1*U2*M2/M6 
1810  T7  =  I6*H*U2*M2/M6 
1820  T4  =  T3*U2*M2/M6 
1830  U3  =  U2*M2/W1 
1840  IF  U3>0 THEN  1860 
1850  GOTO  1990 
1860  R3  =  (T4/W1)/(U3‘  .5) 

1870  M  =  J7+J8  +  J6*G 

1880  PI  =  M/W 

1890  M8  =  W1*(P1-P2)  +  T4*P3 

1900  IF  ABS(M8/Wl)<.  001  THEN  2150 

1910  IF  M8>0  THEN  2030 

1920  A  =  A*(l  +  .1/(1  +•  K8))*5 

1930  K7=  K7+  1 

1940  IF  K7  >  40  THEN  2505 

1950  IF  K7<  15 THEN  1970 

1960  K7=  1.  5 

1970  IF  PI  >0  THEN  2000 

1990  GOTO  1522 

2000  K9=  K9+  1 

2010  K5=  K5  +  1 

2020  GO  TO  570 

2030  A  =  A*(l-.1/(1  +  K7))*  .5 

2031  IF  A>  l  THEN  2040 

2032  PRINT  "YOU  WILL  HAVE  TO  DECREASE  SLOPE  FOR  A  SOLUTION” 

2033  A  =  2 

2034  PRINT  "CURRENT  SLOPE  =”;  M5 

2035  GO  TO  2506 
2040  K8  =  K8  +  1 

2050  IF  K8>30 THEN  2505 

2060  IF  K8<  15  THEN  2080 

2070  K8=1.5 

2080  IF  PI  >0  THEN  2110 

2100  GOTO  1522 

2110  K9  =  K9  +  1 

2120  K5=K5  +  1 

2130  IF  K9> 30  THEN  2505 

2140  GO  TO  570 

2150  PRINT  "SLOPE  =  M5;  "PIVOT  =  ”;P2;",”;P3 
2160  PRINT 

2170  PRINT  "CENT.  PRESS  =  ”;  PI 

2180  PRINT  "PURE  HYDRO  LOAD  =  ”;  W3/B 

2190  PRINT  "MIXED  HYDRO  LOAD  =  ”;  W2/B;  "PSI” 

2200  PRINT  "CONTACT  LOAD  =  W4/B; "  PSI” 

2210  PRINT  "TOT.  LOAD  =  ”;W1/B,  "  PSI” 

2220  PRINT  "PURE  HYDRO  FRICT  =  ”;T5/B.  "PSI” 

2230  PRINT  "MIXED  HYDRO  FRICT  =  ”;  T6/B;  "  PSI” 

2240  PRINT  "CONTACT  FRICT.  =  ”;  T7/B; "  PSI” 

2250  PRINT  "TOT.  FRICT  =  ”;  T4/B;  "PSI" 

2260  PRINT  "U-ETA/W  =  U3 

2270  PRINT  "COEF.  OF  FRICTION  =  T4AV1  +  M5 
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2280  PRINT  "SHOE  SPEED  =  ”;U2;"IN/S” 

2290  PRINT  "PRESS.  MOM.  =  ”;  W1*B*(P1-P2)/12;  "FT  LB/IN” 

2300  PRINT  "SHEAR  MOM.  =  T4*P3*B /12;  "FT-LB/IN” 

2310  PRINT  "H1/H2  RATIO  =  ”;  A 

2320  PRINT  "H2  =  ”;H2 

2330  PRINT  "C/H2”;  Cl 

2340  PRINT  "ABRAMOVITZ’S  PHI  =  F 

2350  PRINT 

2360  PRINT  "NO.  ITERATIONS = K9 

2361  IF  N8 <  .5  THEN  2370 

2362  PRINT  "BEARING  MAY  CAVITATE” 

2370  PRINT 

2380  IF  X4  =  0  THEN  2410 

2390  PRINT  "MIXED  LUBRICATION  BETWEEN”;  X4*100;  "AND”;  (X5  + .  00000*100; 

2400  GO  TO  2420 

2410  PRINT  "NO  MIXED  LUBRICATION” 

2420  PRINT 

2505  PRINT  "FOR  HIGHER  LOAD,  DECREASE  SLOPE,  AND  VICE-VERSA” 

2506  PRINT  "NEW  SLOPE?  (INPUT  ZERO  TO  END  THIS  SEQUENCE)” 

2507  INPUT  M5 

2508  IF  M5C1E-6  THEN  2560 
2510  GO  TO  520 

2560  PRINT  "TOTAL  ITERATIONS  THIS  RUN  =  ”;  K5 

2570  PRINT  "CURRENT  VALUES:  RPM  =  ”;  Ul;  "DIA = ”;  Kl;  "IN. ,  VISCOSITY  =  ”;  M2;  "REYNS, 
2580  PRINT  "LENGTH  =  ”;  B;  "IN.” ;  "ASP  HEIGHT =”;  C2;  "IN. ,  YIELD  STRESS =”;  P8 
2590  PRINT  "RATIO,  SHEAR/YIELD =”;  M3;  "CROWN  HEIGHT  =  ”;  H4;  "IN.” 

2600  PRINT  "PIVOT  CENTER = ”;  P2;",  ”;P3 

2610  PRINT  "DO  YOU  WANT  TO  CHANGE  A  PARAMETER  AND  RERUN?” 

2620  INUT  E$ 

2630  IFSEG*  (E$,  1, 1)  =  "Y”THEN  2670 
2640  IFSEG$(E$,  1, 1)  =  "N”THEN  3370 
2650  PRINT  "INCORRECT  RESPONSE,  ANSWER  YES  OR  NO” 

2660  GO  TO  2620 

2670  PRINT  "NEW  PARAMETER?” 

2680  INPUT  D$ 

2690  IF  SEG$  (D$,  1, 2)  =  "RP”THEN  2820 
2700  IF  SEG$  (D$,  1, 1)  =  "D”  THEN  2850 
2710  IF  SEG$  (D$,  1, 1)  =  "V”  THEN  2880 
2720  IF  SEG$  (D$,  l,  1)  =  "L”THEN  2920 
2730  IF  SEG$  (D$,  1, 1)  =  "A”  THEN  2950 
2740  IF  SEG$  (D$,  1, 1)  =  "Y”THEN  2980 
2750  IF  SEG$  (D$,  l ,  2) = "RA”  THEN  3020 
2760  IF  SEG$  (D$,  1, 1)  =  "C”THEN  3050 
2770  IF  SEG$  (D$,  1, 1)  =  "P”THEN  3080 
2780  IF  SEG$  (D$,  1,2)  =  "RU”  THEN  410 

2790  PRINT  "OPTIONS  ARE  RPM,  DIAMETER,  VISCOSITY,  LENGTH,  ASPERITY,” 

2800  PRINT  "YIELD,  RATIO,  CROWN,  PIVOT,  RUN” 

2810  GO  TO  2670 
2820  PRINT  "NEW  RPM?” 

2830  INPUT  U  l 
2840  GO  TO  2670 
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2850  PRINT  "NEW  DIAMETER,  INCHES?” 

2860  INPUT  K1 
2870  GO  TO  2670 

2880  PRINT  "NEW  VISCOSITY,  REYNS  X  IE6?” 

2890  INPUT  M2 
2900  M2  =  M2*lE-6 
2910  GO  TO  2670 

2920  PRINT  "NEW  LENGTH,  INCHES?” 

2930  INPUT  B 
2940  GO  TO  2670 

2950  PRINT  "NEW  ASPERITY  HEIGHT,  IN  ?” 

2960  INPUT  C2 
2970  GO  TO  2670 

2980  PRINT  "NEW  YIELD  STRESS,  PSI  X  IE-3?” 

2990  INPUT  P8 
3000  P8  =  P8*1E3 
3010  GO  TO  2670 

3020  PRINT  "NEW  RATIO,  SHEAR/YEILD?” 

3030  INPUT  M3 
3040  GO  TO  2670 

3050  PRINT  "NEW  CROWN  HEIGHT,  INCHES?” 

3060  INPUT  H4 
3070  GO  TO  2670 

3080  PRINT  "NEW  PIVOT  CENTER,  X,  Y,  FRACTION  OF  SHOE  LENGTH?” 

3090  INPUT  P2,P3 
3100  GO  TO  2670 
3120  PRINT 

3130  PRINT  "THIS  PROGAM  WILL  CALCULATE  THE  PERFORMANCE  CHARACTERISTICS” 
3140  PRINT  "FOR  A  CROWNED  SLIDER  BEARING  OPERATING  IN  THE  MIXED” 

3150  PRINT  "LUBRICATION  REGIME  THE  INPUTS  REQUIRED  ARE  THE  ROTATIONAL” 
3160  PRINT  "SPEED,  THE  SHOE  LENGTH,  THE  VISOCITY  OF  THE  LUBRICANT,” 

3170  PRINT  "THE  MEAN  BEARING  SURFACE  DIAMETER,  THE  ROUGHNESS  ASPERITY” 
3180  PRINT  "HEIGHT,  THE  BEARING  YIELD  STRESS,  THE  RATIO  OF  THE  SHEAR” 

3190  PRINT  "STRESS  TO  YIELD  STRESS,  THE  CROWN  HEIGHT,  AND  THE  PIVOT” 

3200  PRINT  "COORDINATES,  MEASURED  FROM  THE  LEADING  EDGE  OF  THE  SHOE,” 

3210  PRINT  "AND  EXPRESSED  IN  FRACTION  OF  SHOE  LENGTH  ALL  OTHER” 

3220  PRINT  "DIMENSIONS  ARE  INDICATED  IN  APPROPRIATE  INPUT  STATEMENTS" 

3230  PRINT 

3240  PRINT  "AT  THE  END  OF  EACH  RUN,  ANY  OF  THE  ABOVE  PARAMETERS  MAY  BE  ” 
3250  PRINT  "CHANGE  FOR  A  RERUN.  DURING  EACH  RUN,  THE  OPERATOR” 

3260  PRINT  "WILL  OBTAIN  A  TABULATION  OF  THE  SPEED-LOAD  PARAMETER” 

3270  PRINT  "(SHOE  SPEED  X  VISCOSITYABEARING  PRESSURE  X  LENGTH)” 

3280  PRINT  "VERSUS  COEFFICIENT  OF  FRICTION  AND  THE  RATIO  OF  LEADING” 

3290  PRINT  "AND  TRAILING  EDGE  HEIGHTS.” 

3300  PRINT 

3310  PRINT  "THE  OPERATOR  WILL  BE  ADVISED  IF  HE  HAS  SELECTED  AN” 

3320  PRINT  "UNWORKABLE  COMBINATION  OF  PARAMETERS.  NORMALLY,  CROWN” 
3330  PRINT  "HEIGHT/SHOE  LENGTH  RATIOS  GREATER  THAN  0  0001  WILL  LEAD” 

3340  PRINT  "TO  BEARING  INTERFERENCE,  DEPENDING  ON  THE  PIVOT  POINT  ” 

3350  PRINT 
3360  GO  TO  210 


APPENDIX  B 

INTEGRALS  FOR  CIRCUMFERENTIAL  AND  AXIAL  COMPUTATION 


Integral  formula  is  to  be  developed  for 
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First,  special  characteristics  of  the  second  order  polynomial  will  be  considered 
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(ii)  Degenerate  Linear  Function  b2  =  0,  b^  *  0 
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In  the  event  that  bj  -  ^0^2  <  0>  instead  of  Equation  (B-5),  write 

y2  =  -a2  =  2 b0b2  -  b2 


Then  Equation  ib-8)  can  be  rewritten  as 
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remains  real  and  should  be  used  in  computing  Equations  (C-9),  (C-10),  (C-ll)  and 
C-12 
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O  V'l  * 


DTNSRDC  ISSUES  THREE  TYPES  OF  REPORTS 

1.  DTNSRDC  REPORTS.  A  FORMAL  SERIES.  CONTAIN  INFORMATION  OF  PERMANENT  TECH 
NICAL  VALUE.  THEY  CARRY  A  CONSECUTIVE  NUMERICAL  IDENTIFICATION  REGARDLESS  OF 
THEIR  CLASSIFICATION  OR  THE  ORIGINATING  DEPARTMENT. 

2.  DEPARTMENTAL  REPORTS,  A  SEMIFORMAL  SERIES,  CONTAIN  INFORMATION  OF  A  PRELIM 
INARY.  TEMPORARY.  OR  PROPRIETARY  NATURE  OR  OF  LIMITED  INTEREST  OR  SIGNIFICANCE. 
THEY  CARRY  A  DEPARTMENTAL  ALPHANUMERICAL  IDENTIFICATION. 

3.  TECHNICAL  MEMORANDA.  AN  INFORMAL  SERIES.  CONTAIN  TECHNICAL  DOCUMENTATION 
OF  LIMITED  USE  AND  INTEREST.  THEY  ARE  PRIMARILY  WORKING  PAPERS  INTENDED  FOR  IN 
TERNAL  USE.  THEY  CARRY  AN  IDENTIFYING  NUMBER  WHICH  INDICATES  THEIR  TYPE  AND  THE 
NUMERICAL  CODE  OF  THE  ORIGINATING  DEPARTMENT.  ANY  DISTRIBUTION  OUTSIDE  DTNSRDC 
MUST  BE  APPROVED  BY  THE  HEAD  OF  THE  ORIGINATING  DEPARTMENT  ON  A  CASE  BY  CASE 
BASIS. 


